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We present a multimode laser-linewidth theory for arbitrary cavity structures and geometries that 
contains nearly all previously known effects and also finds new nonlinear and multimode corrections, 
e.g. a correction to the a factor due to openness of the cavity and a multimode Schawlow-Townes 
relation (each linewidth is proportional to a sum of inverse powers of all lasing modes). Our theory 
produces a quantitatively accurate formula for the linewidth, with no free parameters, including the 
full spatial degrees of freedom of the system. Starting with the Maxwell-Bloch equations, we handle 
quantum and thermal noise by introducing random currents whose correlations are given by the 
fluctuation-dissipation theorem. We derive coupled-mode equations for the lasing-mode amplitudes 
and obtain a formula for the linewidths in terms of simple integrals over the steady-state lasing 
modes. 


I. INTRODUCTION 

The fundamental limit on the linewidth of a laser is 
a foundational question in laser theory [Mi It arises 
from quantum and thermal fluctuations [a, l7|, and de¬ 
pends on many parameters of the laser (materials, geom¬ 
etry, losses, pumping, etc.); it remains an open problem 
to obtain a fully general linewidth theory. In this pa¬ 
per, we present a multimode laser-linewidth theory for 
arbitrary cavity structures and geometries that contains 
nearly all previously known effects [sj-fl^ and also finds 
new nonlinear and multimode corrections. The theory 
is quantitative and makes no significant approximations; 
it simplifies, in the appropriate limits, to the Schawlow- 
Townes formula © with the well-known corrections. It 
also demonstrates the interconnected behavior of these 
corrections MM, which are usually treated as inde¬ 
pendent. Most previous laser-linewidth theories have 
employed simple models for calculating the lasing modes 
(e.g., making the paraxial approximation). Such simplifi¬ 
cations, though appropriate for many macroscopic lasers, 
are inadequate for describing complex microcavity lasers 
such as 3d nanophotonic structures or random lasers with 
inhomogeneities on the wavelength scale [T5l - [l8l| . We 
base our theory on the recent steady-state ab-initio laser 
theory (SALT) [Uli^, which allows us to efficiently solve 
the semi-classical laser equations in the absence of noise 
for arbitrary structures [^. We treat the noise as a 
small perturbation to the SALT solutions, allowing us 
to obtain the linewidths analytically in terms of simple 
integrals over the steady-state lasing modes. Our SALT- 
based theory is ab initio in the sense that it produces 
quantitatively accurate formulas for the linewidths, with 
no free parameters, including the full spatial degrees of 
freedom of the system. Hence, we will refer to this ap¬ 
proach as the noisy steady-state ah-initio laser theory 
(N-SALT). 


Our derivation (Secs. III-V) begins with the Maxwell- 
Bloch equations (details in appendix A), which cou¬ 
ple the full-vector Maxwell equations to an atomic 
gain medium (^ . combined with random currents (in 
Sec.|TV| whose statistics are described by the fluctuation- 
dissipation theorem (FDT) [231 - 1?^ . In the presence 
of these random currents, the amplitudes of the las¬ 
ing modes evolve according to a set of coupled ordi¬ 
nary differential equation s ( ODEs), which have been 
called “oscillator models” 1^1^ or “temporal coupled¬ 
mode theory” (TCMT) [soHl in similar contexts. In 
their most general form, our N-SALT TCMT equations 
(Sec. mil) have the form of oscillator equations with 
a non-instantaneous nonlinear term that stabilizes the 
mode amplitudes around their steady-state values. The 
non-instantaneous nonlinearity arises since the atomic 
populations respond with a time delay to field fluctua¬ 
tions; this corresponds to the typical case of “class B” 
lasers in which the population dynamics cannot 

be adiabatically eliminated. We are able to show analyti¬ 
cally that the resulting linewidths of the lasing peaks are 
identical to the results one obtains for a simplified model 
with instantaneous nonlinearity [2^ [2^ , which describes 
the (less common) case of “class A” lasers, in which the 
population dynamics are adiabatically eliminated. As 
expected, however, in certain parameter regimes the full 
non-instantaneous model can exhibit side peaks along¬ 
side the main lasing peaks [1^, arising from relaxation 
oscillations(Sec. IVlC). 

By solving the N-SALT TCMT equations, we obtain a 
simple closed-form matrix expression for the linewidths 
and multimode phase correlations (Sec. |V|, generaliz¬ 
ing earlier two-mode results that used phenomenological 
models . This gives a multimode “Schawlow-Townes” 
relation (Sec. IVII C). where the linewidth of each lasing 
mode is proportional to a sum of inverse output powers 
of the neighboring lasing modes. The theory is valid well 
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random currents 



FIG. 1. (Color online) Schematics illustrating linewidth 
physics, (a) Photonic-crystal (PhC) laser cavity emit¬ 
ting radiation from the lasing mode at frequency cjq, per¬ 
turbed by random currents, (b) The squared amplitude is 
stabilized around a^. Below (above) Oq, the medium provides 
light amplification (attenuation), (c) Phasor diagram for the 
complex field amplitude: a circular oscillation (with |a| = ao) 
for the noise-free mode and a perturbed path for noise-driven 
mode. Noise drives small amplitude fluctuations and possi¬ 
bly large phase drifts, (d) The lineshape is a Lorentzian ~ 
r/[(ai — cJo)^ + (^r)^], centered around uo with width F. 


above threshold, and whenever a new mode turns on, this 
inverse-power relation produces a divergence due to the 
failure of the linearization approximation near threshold. 
However, we show that this divergence is spurious and 
can be avoided by solviiw the nonlinear N-SALT TCMT 
equations numerically |4fll |. (Our formalism can be ex¬ 
tended to treat the near-threshold regime analytically by 
including noise from sub-threshold modes, as discussed 
in Sec. IVIIB and in Sec. IVIIll ) Sec. IVIlIVlII also present 
several other model calculations that illustrate the dif¬ 
ferences between N-SALT and previous linewidth theo¬ 
ries. Finally, in Sec. IVIIIi we discuss some potential ad¬ 
ditional corrections that will be addressed in future work. 
In a second manuscript [4l| , we also compare the theory 
against full time-dependent integration of the stochastic 
Maxwell-Bloch equations and find excellent quantitative 
agreement with the major results presented here. 

Laser dynamics are surveyed in many sources CHI, 
but it is useful to review here a simple physical picture 
of linewidth physics. A resonant cavity [e.g., light bounc¬ 
ing between two mirrors or a photonic-crystal (PhC) mi¬ 
crocavity as in Fig. [TJa)] traps light for a long time in 
some volume, and lasing occurs when a gain medium is 
“pumped” to a population “inversion” of excited states 
to the point (threshold) where gain balances loss. [Of 
course, this simple picture is modified once additional 
modes reach threshold, or for lasers (such as random 
lasers [13, SI]) in which the passive cavity possesses no 
strong resonances; all of these complexities are handled 
by SALT [H [13 and hence are incorporated into our 
approach.] For simplicity, consider here a laser operating 
in the single-mode regime. Above threshold, the gain de¬ 
pends nonlinearly on the mode intensity jap, as sketched 


in Fig. jljb): increasing the field intensity decreases the 
gain due to depletion of the excited states until it reaches 
a stable steady-state value Oq. (This gain-saturation ef¬ 
fect is called “spatial hole-burning” [J since it can be 
spatially inhomogeneous.) In the absence of noise, this 
results in a stable sinusoidal oscillation with an infinites¬ 
imal linewidth, but the presence of noise, which can be 
modeled by random current fluctuations J [13, IH, , 
perturbs the mode as depicted in Fig. [TJc), resulting in 
a finite linewidth. There are various sources of noise in 
real lasers, but spontaneous emission sets a fundamental 
lower limit on the linewidth Q ; here we will include only 
spontaneous emission and thermal noise. In particular, 
although the squared amplitude is stabilized around Oq 
by the nonlinear gain, the phase (j) of the mode drifts 
according to a random walk (a Brownian/Wiener phase) 
with variance ((/^) ~ F t, and the Fourier transform of 
a Wiener phase yields a Lorentzian lineshape [Fig. [TJd)] 
with full width at half maximum (FWHM) F [2g. The 
goal of linewidth theory is to derive F, ideally given only 
the thermodynamic FDT description of the current fluc¬ 
tuations and the Maxwell-Bloch physics of the laser cav¬ 
ity. 

The most basic approximation for the linewidth (suf¬ 
ficiently far above threshold), usua^ referred to as the 
Schawlow-Townes (ST) formula [1, Q, takes the form 


2P 


( 1 ) 


where P is the output power of the laser, 70 is the 
passive cavity resonance width, and luq is the laser fre¬ 
quency, often approximated to be equal to the real part 
of the passive-cavity resonance pole at w* = wq — * 70 / 2 . 
(A slightly more accurate approximation for the laser 
frequency takes into account the small line-pulling of 
the laser frequency towards the atomic transition fre¬ 
quency [ 13 .) The inverse-power dependence causes the 
famous line-narrowing of a laser above threshold. 

Over the decades, a number of now-standard correc¬ 
tions to this formula were found M, leading to the 
modified ST formula: 



( 2 ) 


First, the gain medium can be thought of, in many re¬ 
spects, as a system at negative temperature T [i^], with 
the limit of complete inversion of the two lasing levels 
corresponding to T —>■ 0“. When only partial inver¬ 
sion is present, the linewidth is enhanced by a factor 
of Usp = [ttL 113 , where N 2 and Ni are the spa¬ 

tially averaged populations in the upper and lower states 
of the lasing transition. We refer to this correction as the 
incomplete-inversion factor (also known as “the sponta¬ 
neous emission factor”). Second, due to the openness 
of the laser system, the modes are not power-orthogonal 
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and the noise power which goes into each lasing mode is 
enhanced ; this correction is known as the Petermann 
factor, and it becomes significant in low-Q laser systems, 
where it is not a good approximation to treat the lasing 
mode Ec as purely real. {Q = is a dimensionless 

passive-cavity lifetime defined in units of the optical pe¬ 
riod [s^.) Note that Ec is the passive-cavity mode [in 
contrast to SALT solutions, which are the modes of the 
full non-linear equations, introduced in (l6|)]. f^dx de¬ 
notes integration over the cavity region. Third, for low-Q 
laser cavities, it is possible that the gain linewidth can 
be on the order of or smaller than the passive cavity res¬ 
onance width 7 o, causing significant dispersion effects as 
the gain is increased to threshold This correction is 
commonly called the “bad-cavity” factor Unlike 

the other corrections mentioned above, the bad-cavity 
factor decreases the laser linewidth. However, very few 
lasers systems are in the parameter regime where this ef¬ 
fect is significant [5l|. Finally, amplitude fluctuations in 
the laser field couple to the phase dynamics, leading to 
a correction known as the “a factor”. For atomic gain 
media, this effect was identified by Lax Q in the 1960’s, 
and for this case it is typically a small correction. For 
bulk semiconductor gain media the effect is large, and 
typically dominates the broadening due to direct phase 
fluctuations ■in this context it is known as the 

“Henry a factor” 111. 

Previous linewidth derivations have taken a number 
of different approaches, making severe approximations 
compared to the solution of the full three-dimensional 
space-dependent Maxwell-Bloch equations in the pres¬ 
ence of noise. Generally speaking, linewidth theories 
can be classified into two categories. The first class in¬ 
cludes methods which solve Maxwell’s equations with a 
phenomenological model for the gain medium and ac¬ 
count for noise spatial and spectral correlations by us¬ 
ing the FDT [I 3 , m, . Typically, these methods do 
not handle nonlinear spatial hole-burning above thresh¬ 
old or multimode effects. These methods, commonly 
used in the semiconductor laser literature, resulted in 
linewidth formulas which included the Petermann (^ . 
bad-cavity incomplete-inversion [ 2 ^, and a fac¬ 

tors M- Most notably, an early work by Arnaud 
derived a single-mode linewidth formula without mak¬ 
ing any simplifying assumptions about the field patterns, 
handling anisotropic, inhomogeneous, and dispersive me¬ 
dia. However, this theory was only applied to very sim¬ 
ple, effectively one-dimensional, homogeneous systems, 
and it was missing hole-burning effects and the a factor. 

The second class of linewidth theories consists of 
scattering-matrix methods ( 1 ^ . [13 . [ 5 ^ . [13, which can 
treat arbitrary geometries without phenomenological pa¬ 
rameters and take into account the effects of spatial hole¬ 
burning. S-matrix theories only have access to the input 
and output fields and, therefore, can only treat the noise 
in a spatially averaged manner and are not able to ob¬ 
tain the a factor rigorously. However, they obtain all of 
the other corrections to the single-mode linewidth. In 


particular, the recent S-matrix approach by Chong et 
al. [l3 . [l3 | takes advantage, as we do, of the ab-initio 
computational approach of SALT, and hence has the 
potential to treat arbitrary geometries and spatial hole¬ 
burning effects. (We reduce our results to the most re¬ 
cent scattering-matrix linewidth formula 0 in appendix 
D.) Note that in practice, S-matrix methods require a 
substantial independent calculation beyond SALT to ex¬ 
tract the linewidths, whereas our approach obtains the 
linewidths immediately from SALT calculations (or any 
other method to obtain the steady-state lasing modes) 
by simple integrals over the fields. 


Our derivation of N-SALT, being based on the SALT 
solutions, has a similar regime of validity. For single¬ 
mode lasing^ SALT and N-SALT are essentially exact, 
relying only on the rotating-wave approximation and on 
the laser being sufficiently far above threshold. For mul¬ 
timode lasing, those theories require two additional dy¬ 
namical constraints [l3, 1131 : the rates associated with 
population dynamics must be small compared to both the 
dephasing rate of the polarization and the lasing mode 
spacing (roughly, the free spectral range). The former 
constraint is satisfied in all solid-state lasers, whereas the 
latter requires a sufficiently small laser cavity. The ac¬ 
tual size depends both on details of the cavity and of 
the gain medium used, but the appropriate limit is re¬ 
alized in many complex lasers of interest. When these 
frequency scales are not well-separated, the level popu¬ 
lations are not quasi-stationary, and multimode SALT 
will initially lose accuracy and eventually fail completely 
(since multimode lasing becomes unstable [13 )• More¬ 
over, while the average (SALT) behavior is unaffected 
by non-lasing poles, they do affect the noise properties, 
and N-SALT in its current form only accounts for a finite 
number of poles in the Green’s function (appendix A.2). 
[We only include lasing poles (i.e., poles on the real axis), 
but extension to include non-lasing poles, which deter¬ 
mine the amplified spontaneous emission (ASE) [13 [13 , 
will be straightforward ISec. IVIlH) ]. As noted above, the 
linewidth formula additionally assumes that the laser 
is operating far enough above threshold that amplitude 
fluctuations are small compared to the steady state am¬ 
plitudes (i.e., \a{t)\ ~ oo in the notation of Sec. IV|) . 
Hence, our formula does not describe the linewidth near 
the lasing thresholds. Our perturbation approach takes 
into account only the lowest-order correction to the com¬ 
plex modal amplitude a(t) and neglects higher-order cor¬ 
rections to the frequency loq and spatial pattern Eo(x) 
[see Eq. (O]. Moreover, we neg lect non-Lorentzian cor¬ 
rections to the lineshape [6fl-[6J| fSec. IIVI) . In the follow¬ 
ing section we present our generalized linewidth formula 
in the single-mode regime ([3]) and compare it with tradi¬ 
tional linewidth theories. 
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II. THE N-SALT LINEWIDTH EORMULA 


Symbol Traditional Generalized 


Our main result is a multimode linewidth formula 
which generalizes ([2]). In the multimode case, the re¬ 
sult takes the form of a covariance matrix for the phases 
of the various modes, which is presented in (I36I37I) of 
Sec. |Vl In the single-mode case, the N-SALT linewidth 
formula takes the simple form: 

r=^^-n.p-A-B-(l + 52). (3) 


The modified correction factors (marked by tildes) are 
defined in Table. HI As can be seen from the table, those 
factors generalize the traditional expressions by taking 
into account both spatial inhomogeneity and nonlinear¬ 
ity. Since the generalized factors depend on the SALT 
permittivity e, mode profile Eo(x), and frequency wq, 
one can no longer regard the effects of cavity-openness, 
nonlinearity, and dispersion as separate multiplicative ef¬ 
fects. In this sense, our formula demonstrates the in¬ 
termingled nature of the linewidth correction factors, as 
previously introduced in MM , but here demonstrated 
in a new level of generality. We denote by J dx inte¬ 
gration over all space, for any number of spatial dimen¬ 
sions. We use the shorthand notation for vector prod¬ 
ucts |Eop = Eg • Eg and Eg = Eg • Eg, where the 
latter unconjugated inner product appears naturally be¬ 
cause of the biorthogonality relation for lossy complex- 
symmetric systems [^. Ime(x) denotes the imagi¬ 
nary part of the nonlinear steady-state permittivity ([5|), 
which is negative/positive in gain/loss regions. The out¬ 
put power P is related to the SALT solutions by in¬ 
voking Poynting’s theorem, which one can use to show 
that P oc Jpdx [— Ime(x)]|Eg(x)p. We use Jpdx to de¬ 
note some volume which contains the gain medium. The 
choice of the volume is somewhat arbitrary; e.g., inte¬ 
grating over the cavity region corresponds to the output 
power at the cavity boundary [^. Note, however, that 
this arbitrariness in the choice of the volume is not a 
general feature of our formula. After substituting the 
relevant expressions from Table. U into ®, the integrals 
which contain fp dx cancel, resulting in an expression for 
the linewidth only in terms of integrals over the entire 
space. The effective inverse temperature /3(x) is deter¬ 
mined by the inhomogeneous steady-state atomic popu¬ 
lations A^i(x) and A2(x), and is defined as 


/3(x) = -—In 


^2(x) 


(4) 


In regions where the gain medium is pumped sufficiently 
to invert the population, /3(x) is negative; in regions 
where the pump is too weak to invert, /3(x) will be posi¬ 
tive [and still given by (|4])]; and in unpumped regions, Eq. 
dHI) will simply reduce to the equilibrium temperature of 
the surrounding environment ik^T) The quantities 
Ni{x) and N 2 {x) are an output of the SALT solution 
in the absence of noise. The spatially dependent expres¬ 
sion inside the square brackets in the definition of ngp 


70 

cavity decay 
rate 


B 

bad cavity 


70 


J dx (a;oIn^£)Eo^ 


dx sEo 


^sp 

incomplete 

inversion 

N2 

J 

icoth(«^)-|]lme|Eo|2 

As - Ai 


j rfa:Iine|Eop 


K 

f dx\Ecf 
Jc 

2 

J dxlm € |Eo|^ 

2 

Petermann 

J dxE'^ 


J dx Im € Eo^ 



7J_ + - 


j dxsEg 


j dx Eq 




ujQ ds 


2 du>Q J 


amplitude-phase 

coupling 


7_L 


Im C 

Re C 


C 

nonlinear 

coupling 




J dx(e+^-§^'^ Eg 


TABLE I. Traditional and new linewidth correction factors 
for the single-mode linewidth formulas (I2l3l) . 


in Table. H] generalizes the spatially averaged incomplete- 
inversion factor • That can be seen by noting that 

5 Coth(^^) — i = {exp[hujl3] — 1)“^ = riB, wher e nn 
is the usual Bose-Einstein distribution function . 

(For gain media, it is sometimes convenient to introduce 
the positive spontaneous-emission factor rigp = —ub [Z^- 
Note that this definition ensures that the generalized 
incomplete-inversion factor is always positive.) The ^ 
factor subtracted from the hyperbolic cotangent was dis¬ 
cussed in and we give a simple classical explana¬ 
tion for it in appendix E. If standard absorbing layers 
are used to implement outgoing boundary conditions in 
the SALT solver and the temperature of the ambient 
medium is assigned to these layers, then the N-SALT for¬ 
mula includes the effect of incoming thermal radiation. 
A generalized Petermann factor which formally resem¬ 
bles K appeared in previous work by Schomerus (in 
his expression for the Petermann factor for TM modes 
in two-dimensional dielectric resonators). However, the 
earlier formula is expressed in terms of passive resonance 
scalar fields, whereas our correction contains 3d nonlin¬ 
ear SALT solutions. Finally, 5 is a generalized a factor, 
defined explicitly in Sec. |V] ([30]). For atomic gain media, 
the traditional factor is expressed in terms of the atomic 
transition frequency Wa and decay rate of the atomic po¬ 
larization 7 j_. In the current work we will only evalu¬ 
ate the atomic case, although the general expression in 
terms of the non-linear coupling C should also apply to 
the semiconductor case. 
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The N-SALT formula (jS)) reduces to the traditional 
formula ([2]) in some limiting cases. Let us consider, for 
simplicity, a Id Fabry-Perot laser cavity of length L sur¬ 
rounded by air (i.e., Ime = 0 outside the cavity region). 
Let us assume also that the laser is operating not too 
far above the threshold and is uniformly pumped, hence 
Im e and P are nearly constant inside the cavity. In this 
limit, all the integrals in Table. U can be approximated 
by reducing the integration limits to the cavity region; 
terms which contain integration over the imaginary part 
of the permittivity are non-zero only within the cav¬ 
ity region (e.g., / c?xIm£r|Eop becomes Ime j^dx |Eop); 
while terms of the form J da; e Eq can be written as the 
sum of the cavity contribution e dx Eq and the sur¬ 
rounding medium contribution dx Eq, where the lat¬ 
ter is negligible for Lwq ^ I, as shown in appendix D and 
in 0 (here and throughout the paper, we are setting 
c = 1). Using this approximation, it is immediately ap¬ 
parent from Table. U that the incomplete-inversion factor 
reduces to the traditional expression. The generalized 
Petermann factor reduces to the traditional factor in the 
limit of a high-Q cavity, where the threshold lasing state 
E(x) is approximately the same as the passive resonance 
state Ec(x). In order to simplify the remaining terms, re¬ 
call that the lasing threshold is reached when gain in the 
system compensates for the loss. For weak losses (small 
Ime/e) that can be treated by perturbation theory, the 

threshold condition is 70 = 0 and, therefore, the 

generalized decay rate reduces to 70 (one can thereby see 
that the Schawlow-Townes formula neglects nonlin¬ 
ear corrections to 70 , as was also shown in [l^)- Next, let 
us discuss the generalized bad-cavity factor, which sim¬ 
plifies to ^1 -I- after reducing the integration 

limits. In order to show that it agrees with the tradi¬ 
tional factor, we need to show that ~ 2 ^- 

steady-state effective permittivity, as used in SALT the¬ 
ory (appendix A.I), is 


e(x) = ec(x) 


7±^(x) 

Wo - Wa -I- i7_L ’ 


(5) 


where £c is the passive permittivity and the second 
term is the active nonlinear permittivity due to the gain 
medium. The population inversion D(x) = Af2(x) — 
A^i (x) is generally spatially varying above threshold due 
to spatial hole-burning. Since we assume here that we are 
close to threshold and that the pumping is uniform, the 
inversion is also uniform in space and near its threshold 
value. If one assumes, additionally, that the detuning 
of the lasing frequency from atomic resonance is small 
(|wo - Wa| < 7 -l), one obtains Ri Finally, we 

show in Sec. IVII A that our a reduces to the known ao in 
homogeneous low-loss cavities, so that all factors of the 
corrected ST formula are recovered in this limit. (Note 
that line-pulling effects which may modify the lasing fre¬ 
quency Wq are handled by SALT.) 

In the next section, we present the TCMT equations 
which are used in this paper to derive the N-SALT 


linewidth formula ®, but which may also be used to 
extract more information on laser dynamics away from 
steady state. 


III. THE N-SALT TCMT EQUATIONS 

In the absence of noise, the electric field of a laser op¬ 
erating in the multimode regime is given by the real part 
of Eo(x,t), where 

Eo(x,t) = ^E^(x)a^oe“*‘^'‘*, (6) 

and the laser has zero linewidth. (This assumes, of 
course, that there exists a steady-state multimode solu¬ 
tion of the nonlinear semi-classical lasing equations fiol . 
[ 2 ^.) The modes E^(x) and frequencies can be 
calculated using SALT, which solves the semi-classical 
Maxwell-Bloch equations in the absence of noise. (SALT 
has been generalized to include multi-level atoms [t^, 
multiple lasing transitions, and gain diffusion [^; any of 
these cases can thus be treated by N-SALT with minor 
modifications, but we focus on the two-level case here.) 
The linewidth can now be calculated by adding Langevin 
noise, as described below. 

In the presence of a weak noise source, the electric 
field can be written as a superposition of the steady- 
state lasing modes with time-dependent amplitudes (t) 
which fluctuate around 

E(x, t)=Y^ E^(x)a^(t)e"*‘^'‘‘. (7) 

In principle, the sum in © should also include the 
non-lasing modes since the set of lasing modes by it¬ 
self does not form a complete basis for the fields. Non¬ 
lasing modes contribute to amplified spontaneous emis¬ 
sion (ASE), which has a significant effect on the spectrum 
near and below the lasing thresholds [13, and will be 
treated in future work. 

In appendix A, we derive the N-SALT TCMT equa¬ 
tions of motion for a^{t) starting with the full vectorial 
Maxwell-Bloch equations. We show that the noise-driven 
field obeys an effective nonlinear equation which, in the 
frequency domain, takes the form 

[V X V X —a;^e(a;, a)] E(x,a;) = Fs(x,a;), (8) 

where the carets denote Fourier trans¬ 
forms [e.g., E(x,t) = Jq'^ dcLie“®‘^‘E(x,a;)]. Sponta¬ 
neous emission is included via the stochastic noise 
term Fs(x, w) (quantified in Sec. IIVI) . and the effective 
permittivity e(oj,a) (derived in appendix A.2) is given 

by 


e(uj, a)E(x, w) = [ecSm + * “m] E^(x) 


( 9 ) 
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where the asterisk denotes a convolution. The second 
argument of e(w, a) denotes the implicit dependence of 
e on the modal amplitudes through the inversion D. 
The effective permittivity dH) can be decomposed into a 
steady-state-amplitude dispersive term and a nonlinear 
non-dispersive term (similar in spirit to [zl). The key 
point here is that, to lowest order, there are two correc¬ 
tions to the permittivity in the presence of noise: the 
dispersive correction due to any shift in frequency at the 
unperturbed amplitudes a^O) and the nonlinear correc¬ 
tion due to any shift in amplitude at the unperturbed 
frequency. (Shifts in frequency are small because only fre¬ 
quency components within the mode linewidths matter, 
while shifts in amplitude are small because of the stabiliz¬ 
ing effect of gain feedback.) The coupling between these 
two perturbations is higher order and is hence dropped, 
which greatly simplifies the analysis. 

Substituting the permittivity expansion (derived ex¬ 
plicitly in appendix A.3) into Maxwell’s equation ([ 8 |), we 
find that the noise-driven field obeys the linearized equa¬ 
tion 

[V X V X -u}^e{uj,ao)] E(x,w) = Fnl(x,w) + Fs(x,a;), 

( 10 ) 

i.e., the dispersive permittivity which appears on the left- 
hand side of m is evaluated at the steady-state am¬ 
plitude oq. The nonlinear non-dispersive term Fnl [de¬ 
fined explicitly in (IA23I) ]. which corresponds to amplitude 
fluctuations at the unperturbed frequency, appears as a 
restoring force on the right-hand side. The noise-driven 
field F(x, w) is found in appendix A.4 by convolving the 
linearized Green’s function with the source terms Fnl 
and F 5 . Finally, the N-SALT TCMT equations are ob¬ 
tained by transforming the noise-driven field back into 
the time domain. 


A. Time-delayed multimode model 


We find that, in the most general case, the TCMT 
equations take the form 


1 / 

7 (x) 


O-fj. + fii- 

( 11 ) 


Comparing m and (113, one can see that the first term 
on the right-hand side of dm is related to the nonlin¬ 
ear restoring force Fnl, and the Langevin noise /^(t) is 
associated with Fg. 

The nonlinear coupling coefficients c^j/(x) [derived in 
(IA33|) ] correspond to local changes in the nonlinear per¬ 
mittivity with respect to intensity changes in each of the 


modes 


— lUJ 


2 de{u}^ 


F2 


C-LLU - 




( 12 ) 


where we have introduced a shorthand notation for the 
derivative in the denominator (ujf.e)' = . This 

modal coupling in the fluctuation dynamics comes about 
because of saturation of the gain: a fluctuation in mode 
fjL affects the amplitudes of all the other modes v. 

The N-SALT TCMT equations are nonlocal in time 
because the atomic populations are not in general able to 
follow the field fluctuations instantaneously and, instead, 
respond with a time delay determined by the local atomic 
decay rate 7 (x), given by 


■<i^) = 711 (1 + E . 

(13) 

The second term in (jl3l) is precisely the local enhance¬ 
ment of the atomic decay rate due to stimulated emission 
in the presence of the lasing fields. (A simplified spatially 
averaged enhancement of the atomic decay rate was pre¬ 
viously discussed in [t^.I 

The Langevin force is the projection of the spon¬ 
taneously emitted field onto the corresponding mode 
F^ [ 2 ^. Defining F^(t) = Fse®'^'"*, the Langevin force 

/a* is 

i I dxEfj, ■ F^(t) 

Mt) = —r -■ ( 14 ) 

J 

The full N-SALT TCMT equations (HID describe the 
most typical situation in laser dynamics of a “class 
B” laser [35l-l3^. in which the polarization of the gain 
medium can be adiabatically eliminated but the popula¬ 
tion dynamics is relatively slow and cannot be so elimi¬ 
nated. However, much of the basic linewidth physics can 
be extracted from the limit when the population dynam¬ 
ics is also adiabatically eliminable, which describes “class 
A” lasers. Since the mathematical analysis is simpler in 
this limit, we will begin the spectral analysis in Sec. |V] 
with the latter model. We discuss this limit, which we 
refer to as the “instantaneous model,” in the following 
section. 


B. Instantaneous single-mode model 

When the population relaxation rate 7 (x) is (every¬ 
where) large compared to the dynamical scales deter¬ 
mining a^(t), the exponential terms in (HID act like 6 
functions. After the spatial integration, and specializing 
in this section to the single-mode case, we obtain the sim¬ 
ple nonlinear oscillator model driven by a weak Langevin 
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force f{t): 

d = C {al-\a\‘^)a +f, (15) 

where C = J dx c(x) is the integrated nonlinear coupling. 
This instantaneous nonlinear oscillator model was previ¬ 
ously introduced by Lax [§, [^, and has been used ex¬ 
tensively in linewidth theories 0. The N-SALT approach 
enables computing the model’s parameters ab initio, tak¬ 
ing full account of the spatial hole-burning term and the 
vectorial nature of the fields [including multimode ef¬ 
fects, when generalizing m to the multimode regime]. 
Also, our approach shows that this well-known model 
can be explicitly derived from the more general (non- 
instantaneous) model, presented in the previous section. 
Above the lasing threshold, oq > 0 and Re[C] > 0, 
and the system undergoes self-sustained oscillations with 
a stable steady state at jaj = ag, as demonstrated in 
Fig. mb). In fact, near threshold one can show that Re[C] 
is approximately the threshold gain, which balances the 
cavity loss k. Hence the dynamical scale of a(t) is of or¬ 
der K, which must then be much smaller than 7 (x) for 
the instantaneous model to hold; this is the standard dy¬ 
namical condition for class A lasers 

The nonlinear term in m and the multimode coun¬ 
terpart in (fT^ are derived rigorously in appendix A, but 
we can motivate the resulting expressions using simple 
physical arguments. The nonlinear term can be viewed 
as a shift in the oscillation frequency, i.e., —iAco = 
C{aQ — jap). Using first-order perturbation theory [ll, 
the frequency shift due to a change in dielectric permit¬ 
tivity Ae is given by 


Au = —u 


dx Ae Eg 


dx(u:le)'gEil 


(16) 


Plugging in the differential of the permittivity due to 
small changes in the squared mode amplitude, Ae k, 
g^(|ap — Og), we find that the coupling coefficient in 
the instantaneous model is 


9 , de 

J dx^wleYf^Bl 


(17) 


This is the single-mode version of (113 integrated over 
space due to rapid relaxation. As we will see, this simple 
result, combined with the spectrum of the Langevin noise 
(section IV), is all that is needed to derive the single¬ 
mode linewidth formula (|3|) (see Section V), and the 
multimode generalization also follows straightforwardly. 
Hence, after analyzing the noise spectrum, we will first 
derive the linewidth within the instantaneous model be¬ 
fore moving on to the more complicated case of the full 
N-SALT TCMT equations. The latter will show that the 
basic linewidth formula is unchanged from that of the in¬ 
stantaneous model except for the addition of side peaks 
due to the relaxation oscillations present in class B lasers. 


IV. THE AUTOCORRELATION EUNCTION OF 
THE LANGEVIN FORCE 


In this section, we express the autocorrelation function 
of the Langevin force 

= (18) 

in terms of the autocorrelation function of the noise 
source F^. It is well known that quantum and ther¬ 
mal fluctuations can be modeled as zero-mean random 
variables, defined by their correlation functions 1^. 
This Rytov picture [2^ is essentially a consequence of 
the central-limit theorem (CLT) [^, which holds 
since the classical forcing Fs is the sum of a large num¬ 
ber of randomly emitted photons. The autocorrelation 
function of Fs can be found by invoking the fluctuation- 
dissipation theorem (FDT), as explained below. 

The probability distributions of the pumped medium 
and the electromagnetic field obey Boltzmann statis¬ 
tics, with an effective local temperature /? defined in 
terms of the atomic inversion |46l | (see definition in 
Sec. HU- Under the typical conditions of local ther¬ 
mal equilibrium |23l427j . dissipation by optical absorp¬ 
tion must be balanced by spontaneous emission from 
current fluctuations J(x, t). One can then apply the 
FDT for the Fourier-transformed forcing Fs(x, w) = 
—47riwJ(x,w) [SOj j: 

(Fs(x,w)Fs(x',w')^ = 

e{x, uj) coth ^ (5(x — x')d(w — uj'). 

(19) 

Using this result, we calculate the autocorrelation of 
the Langevin force [i.e., the Fourier transform of m, 
defined as /^(w) = ^ /g'^ dt e*“‘/^(t)] and we obtain 

(7^(tu)7:(w')) = Moj)S{io - (20) 


where the frequency-domain autocorrelation coefficient is 




J dx |E^|^Ime(ai) 

— coth 

2 

' hujfd'' 



J 

f dxElixley^ 



( 21 ) 

The I factor subtracted from the hyp erbolic cotangent 
is explained in appendix E and in [T^. 

The time-domain diffusion coefficient can be found 
directly from (1211) taking the inverse Fourier transform. 
For the common case of a small linewidth, Ime(a;) and 

coth ( ) are nearly constant for frequencies within the 


linewidth. [This means, essentially, that the Langevin 
force (t) can be treated as white noise]. Consequently, 
one can approximate the diffusion coefficient in (1211) by 
its value at With this simplification, the time-domain 












diffusion coefficient in (flSl) is conveniently given by = 

More generally, however, including this frequency de¬ 
pendence corresponds to temporally correlated fluctua¬ 
tions, leadi^ to non-Lorentzian corrections to the laser 
lineshape (hOl - l^ . These “memory effects” can be ad¬ 
dressed using our approach (as discussed in Sec. IVIIII) 
and we plan to include them in future work. 


V. THE LASER SPECTRUM 

In this section, we calculate the laser spectrum using 
the N-SALT TCMT equations (I11I15I) and the noise au¬ 
tocorrelation function (I20I21I) . We begin by showing that 
the phase of the lasing mode undergoes simple Brownian 
motion; consequently, the laser spectrum is a Lorentzian, 
with a width given by the phase-diffusion coefficient. 
In Sec. ElA, we calculate the phase-diffusion coefficient 
(hence the linewidth) for the instantaneous model (fT51) 
and in Sec. Eb, we outline the analysis for the time- 
delayed model dm, leaving the details of the derivation 
to appendix B. More accurately, the spectrum of the 
time-delayed model consists of a central Lorentzian peak 
at the lasing resonance frequency and additional side 
peaks due to relaxation oscillations, which are present 
in class B lasers. The latter side peaks are the subject of 
Sec. EC. 


the phase variance increases linearly in time. An oscil¬ 
lator with Brownian phase noise has a Lorentzian spec¬ 
trum [8l|, and one can reproduce that result briefly as 
follows. The laser spectrum 5a(w) is given by the Fourier 
transform of the autocorrelation function of a(t): 


(a(t)a*(0)) « Oo (e 


2 ( = ^2 -i((0(t)-0(o))=') 


/ 


(26) 

For a Wiener phase, whose variance is ((<(>(t) — ())(0))^) = 
r|t|, the Fourier transform of (051) is a Lorentzian whose 
central-peak width is F 2^. In passing from the first 
to second step in (051) . one neglects direct amplitude- 
fluctuation contributions (which are decoupled from the 
phase) as these only introduce broad-spectrum back¬ 
ground noise, but do not affect the linewidth of the laser 
peak (we return to this point in Sec. EC). In passing 
from the second to the third step, one assumes that the 
phase is a Gaussian normal variable, which is justified as 
a consequence of the CLT. 

It is well known that also in the general case of B ^ 0, 
the phase is a Wiener process, with a modified diffusion 
coefficient [HI- In order to calculate the phase variance 
explicitly, we solve (I24I25I) and obtain 




m=j 

ao(/>(t) = —B J 5{t')dt' + J fi{t')dt'. 


(27) 

(28) 


A. Instantaneous single-mode model 

The complex mode amplitude a(t) can be written in 
polar form as 

a(t) = [ao + 5(t)] (22) 

tto is the steady-state amplitude, while i5 and (j) are 
real amplitude and phase fluctuations. Substituting the 
modal expansion ((2^ in (fTSl) . defining 


Substituting (EZl) into (1^ . using the autocorrelation 
function of / (EU, and performing the integration, one 
obtains that the phase variance in the long-time limit 
is (((/)(t) - (/)(0))2) = ^ (1)^^ |t| (where terms 

growing more slowly than |t| were neglected, as explained 
in greater detail in appendix B). Therefore, the linewidth 
is 

r = ^(i + 5^), (29) 


A = 2alRe C 

B = 2allmC, (23) 

and keeping terms to first order in S/qq, we obtain 

(5 = -A5 + /r, (24) 

ao0 = -BS + fi, (25) 

where /r = Re {/} and /i = Im {/}. We check the ap¬ 
proximation of |<5| <C oq a posteriori and we find that 
it generally holds (as was also shown in 0), except 
near threshold {oq 0), which is a case we discuss in 

Sec.EIlC. 

When the nonlinear coupling coefficient is real {B = 0), 
it is evident from (1251) that the phase undergoes simple 
Brownian motion (i.e., it is a Wiener process) and hence 


where we have defined the generalized a factor: 


R _ ImC 
A ~ ReC’ 


(30) 


with the nonlinear coefficient C defined in (EZD. Substi¬ 
tuting the autocorrelation function (|2T]) in (1291) and us¬ 
ing Poynting’s theorem to relate Oq to the output power 

p _ Jp da; (—Ime(x))|Eo(x)p 0, we obtain the 
single-mode linewidth formula ([3]). From (|29)) . it is evi¬ 
dent that the Schawlow-Townes, Petermann, bad-cavity 
and incomplete-inversion factors are all included in the 

term , and generally cannot be separated into the tra- 
2“0 ,—, 

ditional factors of ([2|) |14| . 

When the nonlinear coupling coefficient is complex 
(i.e., when R 0), the resonance peak is not only broad¬ 
ened but is also shifted [HI- The shift in center frequency 
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FIG. 2. (Color online) Simulated spectrum Sa(ca) of the in¬ 
stantaneous model m with ReC = 10, noise autocorrelation 
coefBcient 7? = 0.1 and three values of 5: 10 (blue), 5 (red), 
1 (yellow) {C,R,Sa, and uj are given in arbitrary frequency 
units). The noisy signal is the simulation result and the black 
curves are Lorentzian lineshapes with widths F and center 
frequency shifts given by ((2913111 . 


is found by keeping second-order terms in S/ao and cal¬ 
culating the average phase drift: 


6u} = {(j)) 


RB 

Aa^A' 


(31) 


approximation for a uniformly pumped class B laser op¬ 
erating near threshold). In this case, the single-mode 
time-delayed model takes the form 

a = c(^-foJ a +f, (33) 

where C = J dx c(x) is the integrated nonlinear coupling 
and c(x) is defined in (fT^ . This integro-differential equa¬ 
tion can be turned into a first-order ODE by using the 
modal expansion from Sec.|VlA: a = (oq + d)e®‘^, keeping 
terms to first order in S/ao, and introducing the variable 

f (34) 

Then, (I33I34I) can be recast in the form v = IKv-bf, where 
V = {(5, ao(^,0. 

However, most generally, the spatial dependence of 
7 (x) cannot be neglected. The time-averaged deviation 
^(x, t) is therefore spatially dependent, and one obtains 
an infinite-dimensional problem. To simplify the alge¬ 
bra, we discretize space [e.g., discretizing (fTTl) into a Rie- 
mann sum over sub-volumes 14 ] and recover the con¬ 
tinuum limit at the end. This yields the discrete-space 
multimode model: 


An identical formula was derived in in a phenomeno¬ 
logical instantaneous model. 

Fig.[2]shows the spectrum of the instantaneous model, 
which is obtained by numerically solving (USD using a 
stochastic Euler scheme [s^. Introducing the notation 
R{a) = C{aQ — |ap)a and discretizing time as a{nAt) « 
Gn, the Euler update equation for the n-th step is 

Gn = a„-i + T (a„_i) At -I- Vi? At C, (32) 

where C is a gaussian random variable of mean 0 and 
variance 1, i.e., (/ G N{0,1). [For the data presented 
in Fig. 121 At was decreased until the simulation re¬ 
sults converged. In later sections (Fig. |3D, we imple¬ 
mented a fourth-order Runge-Kutta method in order to 
achieve convergence]. The simulated spectra (noisy col¬ 
orful curves) match the predicted Lorentzian lineshapes 
(solid black curves), which are calculated using (|29I31D . 
As 5 increases, the linewidths are broadened and the 
center frequencies are shifted. 


= 

(ik _ |a,(i')|2)^ a, + f„ 

vk ^ ' 

(35) 

where the discretized nonlinear coupling coefficients are 
= /i 4 dxcy,,,{yi) (so that O'/,,), jk is the 

relaxation rate at the fc’th spatial point and a^o is the 
steady-state amplitude of mode v. 

In appendix B, we study the statistical properties 
of the solutions to (1351) . We introduce the the M- 
dimensional vectors whose entries are = (where 

M is the number of active lasing modes) and we calcu¬ 
late the covariance matrix (<i>^(t)<I’i.(0)). We find that 
the result is independent of the relaxation rates or the 
discretization scheme: 

^(0)) = Q + BA-1 ^ |t|. (36) 


B. Time-delayed multimode model 

We now turn to the laser spectrum produced by the 
time-delayed model, where the nonlinearity is dependent 
on the modal amplitudes at previous times. Although 
we calculate the linewidth of the full time-delayed N- 
SALT TCMT equations (fTTl) in appendix B, we begin 
this section by considering the simplified case of a spa¬ 
tially homogeneous medium 7 (x) Ri 70 (this is a good 


The matrices A and IB correspond to the real and imagi¬ 
nary parts of the coupling matrices, with entries = 
2u^o^i^oF.e[C^i/] and B^u — 2(j^QQ.,yQlm[C^i/]. IR is the 
autocorrelation function of the Langevin force vector f 
[defined in (l2T]l ]. The diagonal of this matrix, divided 
by |t| and by the squared modal amplitude, gives the 
generalized linewidths 


r/.= 


1 





) 



(37) 

















10 



FIG. 3. (Color online) Simnlated spectrum of the time- 
delayed model (I33II with ReG = 10 and R = 0.1 (in arbitrary 
frequency units) at six values of 70 (using a base 10 logarith¬ 
mic scale for the y axis). The noisy signal is the simulation 
result and the black curves are Lorentzian lineshapes with 
widths given by (1291) . 


Therefore, the generalized a factor (which is responsible 
for linewidth enhancement due to coupling of amplitude 
and phase fluctuations) is given by 


Ru 



(38) 


In the single-mode case {M = 1), this matrix formula re¬ 
duces to the single-mode linewidth of the instantaneous- 
model: ^(1 + (§)^) [ (|29I30I) in Sec. ElA]. 

The linewidth in the time-delayed (class B) model is 
precisely the same (neglecting side peaks) as in the in¬ 
stantaneous (class A) model. While this result was de¬ 
rived for single-mode class B semiconductor lasers using a 
phenomenological rate-equation framework [t^ , we prove 
that this is generally the case in the multimode inhomo¬ 
geneous regime. Naively, one might expect to obtain dif¬ 
ferent linewidths due to the longer time over which the 
fluctuations can grow. However, in appendix B we ob¬ 
tain a linewidth expression which is independent of the 
relaxation-oscillation dynamics, which demonstrates that 
there is a cancellation of two competing processes: as yy 
decreases, amplitude fluctuations grow, but they are also 
averaged over longer periods of time so that their effect 
is smaller. 

Fig. [3] presents the simulated spectrum of the time- 
delayed model in the homogeneous -7 limit, which is ob¬ 
tained by numerically integrating (1331) (by applying a 
stochastic Euler scheme, as in Fig. [2]). The width of 
the central peak of the spectrum matches our prediction 
dm, independent of the value of 70 . At intermediate 
relaxation rates, we also observe side peaks in the spec¬ 
trum due to amplitude relaxation oscillations (RO), in 
addition to the central peak. 


C. Side peaks in the time-delayed model 

In class B lasers, amplitude fluctuations relax to steady 
state via relaxation oscillations [ 1 ^ and, consequently, 
give rise to side peaks in the spectrum, in analogy with 
amplitude modulation of harmonic signals. Mathemat¬ 
ically, the oscillation arises from the second-order ODE 


generated by coupling of the 5 and ^ equations (|33I34I) . 
producing the coupled amplitude/gain oscillations. Us¬ 
ing the same methods that we applied to calculate the 
linewidth of the central resonance peak dm, we also 
calculated the full side-peak spectrum in the multimode 
regime. Our formula is derived under the fairly general 
assumption that the central resonance peaks are narrower 
than the side peaks, which is the relevant regime for many 
lasers [t^. Although the derivation uses the same tech¬ 
niques as in appendix B, it is fairly involved and will be 
provided in a subsequent manuscript [s^; we only sum¬ 
marize here. 

As was shown in Sec. IIIIl far above threshold, the 
atomic relaxation rate (1131) is enhanced and can even be 
dominated by the electromagnetic field. This modified 
relaxation rate, and in particular its spatial dependence 
due to hole-burning effects, has important implications 
on the RO spectrum which, to our knowledge, have not 
been treated before. For simplicity, we focus here on the 
case of a = 0. (Note that a factor effects on the RO 
spectrum have been observed and analyzed using a phe¬ 
nomenological homogeneous time-delayed model in |76l|.) 

In order to see how one can obtain a closed-form ex¬ 
pression for the RO spectrum, recall that when calcu¬ 
lating the spectrum of the central resonance peak in 
Sec. El A, we neglected direct amplitude-fluctuation con¬ 
tributions in (EH), i.e., in passing from the first to second 
step, we omitted a term of the form 

((5(f),5(0)) • (39) 


Adding this term in (1331) . one finds that the full spec¬ 
trum consists of an additional term, which is given by 
the convolution of the real-amplitude fluctuation spec¬ 
trum ((5(t)i5(0)) and the spectrum of the central reso¬ 
nance peak. In the instantaneous model, the amplitude 
autocorrelation function {5{t)5[Q)) decays exponentially 
in time [see E3] and the omitted term results in near¬ 
constant background noise. However, in the time-delayed 
model, this neglected term is responsible for the RO side 
peaks. 

For simplicity, consider first a model which can be 
solved straightforwardly; the single-mode homogeneous- 
7 time-delayed model [i.e., 7 (x) Ri 70 and / dxc{'x) = C 
as in (1331) ]. which describes uniformly pumped single¬ 
mode lasers near threshold. Following the discussion 
in Sec. ElB, we can rewrite (I33I34I) as a set of linear 
equations and solve for 5{t), obtaining 



A 


sinh 


-it-t') 


Mt'), 

(40) 


where A = 7 / 7 ^ — 4A 70 . In the limit of well-resolved 
side peaks (e.g., ^ <C 70 “C A), the amplitude autocor- 
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relation function is approximately 


{smo)) 


R 

"2 


sin cos \/^QAt 


7o 


X e 


(41) 

Thus, additional peaks in the spectrum arise at frequen¬ 
cies WRo = Wo ± V^7o with widths 70 . In the high- 
Q limit near threshold, A is proportional to the cavity 
decay rate k, giving the expected behavior for the RO 


frequency. The side-peak amplitudes 
verge in the limit of 70 


70 


di- 


0 (that is, when amplitude 
fluctuations are not small compared to the steady-state 
mode amplitude), but this is also the regime in which our 
analysis of the spectrum (Sec. [V] A-B) breaks down. The 
inset in Fig. shows the simulated spectrum of the ho¬ 
mogeneous time-delayed model (l33|) (the same data was 
also shown in Fig. [31 but we include here the theoretical 
formula for the side-peak spectrum). The exact numeri¬ 
cal solution of (1331) (blue curve) reproduces the analytic 
spectrum prediction of (1421) (red curve). 

In the limits of extremely small/large relaxation rates 
7 o (compared to A), the side peaks disappear. In the 
former limit, they merge with the central resonance peak 
and in the latter case, they merge with the background 
noise. This behavior can be explained by inspection of 
the S and ^ equations (I33I34I) in the appropriate limits. 
When the relaxation rate is very large, the time-delayed 
model reduces to the instantaneous model, which repre¬ 
sents the case where the atomic population follows the 
field adiabatically. In the opposite limit of extremely 
small relaxation, the field follows the atomic population 
adiabatically. In other words, a clear separation of atomic 
and optical time scales will result in the absence of RO 
side peaks. 

In the most general spatially inhomogeneous time- 
delayed model, the full spectrum takes the simple form 


5a(w) = 


-( 5 ) 






1 — dx 


A(x)7(x) 


+ (7 +7(x))2 


h 


A(x)7(x) (7(?i 


-t (7(x)-t 5 )' 


(42) 


where A(x) is the real part of the local nonlinear coupling 
[defined in (1331) ]. 7 (x) is the effective decay rate, and 
F is the central peak linewidth. (This formula is valid 
when the central resonance peak is narrower than the 
side peaks F <C 7 ||.) Like our linewidth formula, this 
formula is easy to evaluate via spatial integrals of the 
SALT solutions. 

While the homogeneous time-delayed model near 
threshold agrees with standard results on relaxation os¬ 
cillations [^, the full model above threshold, combined 
with SALT, is able to include effects not contained in 
other treatments. As the pump is increased far above 
threshold, the effects of stimulated emission strongly in¬ 
crease the atomic relaxation rate, and spatial hole burn- 


(a) 
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FIG. 4. (Color online) Dressed decay rate and the RO spec¬ 
trum based on SALT solutions of a Id PhC laser. Inset: a 
quarter-wave PhC (period a = 1 mm and alternating layers 
with permittivities ei = 16 + 0. li and £2 = 2 + 0. li and thick¬ 
nesses di = and d 2 = a — d\). The center region has 

permittivity ea = 3 + O.li and contains gain atoms with band¬ 
width 7 x = 3 mm~^ and resonance frequency tUa = 25 mm“^. 
(a) Dressed decay 7(x) evaluated using (1131) at five pump val¬ 
ues ( 2 Dth brown, 3Dth blue, 4Dth black, and 5Dth gray), (b) 
Side-peak spectrum Sro(cu) evaluated using (1421) for the hve 
pump values of (a). Inset: full simulated spectrum Ss.{ui) on 
a semi-log scale (of base 10) of the homogeneous time-delayed 
model (1331) with 70 = 0.09, A = 10, B = 0, R = 0.01 (in ar¬ 
bitrary frequency units). The noisy signal is the simulation 
result and the red curve is the theoretical lineshape (HU). 


ing causes that rate 7 (x) to vary substantially in space 
[see mi These two effects cause both a shift and 
a broadening of the side peaks compared to the near¬ 
threshold result. Fig. m shows the dressed decay rate 
7 (x) and the side-peak spectrum S'Ro(a;) [as given by 
the second term of ((131) ]. based on a SALT calculation of 
a one-dimensional photonic crystal (PhC) laser, at four 
different pump values well above threshold. [The pump 
value is controlled via the parameter Dp in (IA3p . and we 
denote the threshold value of Dp by Dth-j This type of 
cavity (depicted in the inset of Fig. |4^) supports a single 
mode at the simulated parameter regime, which is local¬ 
ized near the defect region. (Further discussion of this 
structure is given in Sec. [VllA below.) As can be seen 
from Fig. [4^, the decay rate 7 (x) is enhanced at high 
intensity regions (i.e., near the defect), and it increases 
further as the pump increases. Fig.jljD demonstrates the 
shifting and broadening of the side peaks. 


VI. THE GENERALIZED a FACTOR 


Our TCMT derivation of the linewidth formula yields 
a generalized a factor (l38l) which depends on the eigen- 
modes E^(x) and eigenfrequencies of the full non- 
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FIG. 5. (Color online) (a) The generalized (blue) and traditional (red) a factors of a PhC laser vs. relative detuning Azz = ■ 

Upper inset: quarter-wave PhC geometry (see caption of Fig. SI). The gain parameters are 7 x = 3 mm ^ and a varying cJa- 
Lower inset: intensity distribution of the lasing mode, (b) a factor for an open cavity laser vs. passive permittivity Ec- 
Blue (red): generalized (traditional) a factor. Upper inset: dielectric slab, of permittivity £c, bounded by air on both sides, 
containing gain atoms with Wa = 15 mm“^, yx = 3 mm~^. Lower inset: intensity distribution of the lasing mode. Leftmost 
inset: enlarged segment of the main plot, around £c = 7. 


linear SALT equations. This is an advance over previ¬ 
ous linewidth formulas; the ab-initio scatterin g-m atrix 
linewidth formulas did not obtain an a factor [13l . [a, 
whereas other traditional laser theories that derived a 
factors could not handle the full nonlinear equations . 
Therefore, in the following section, we focus on the gen¬ 
eralized a factor. We compare the generalized and tra¬ 
ditional factors in Sec. EDA, and then we evaluate these 
factors in the single-mode (Sec. IVIIB) and multimode 
(Sec.EDC) regimes. 

A. Comparison with traditional a factor 

Linewidth broadening due to amplitude-phase cou¬ 
pling (that is, the a factor linewidth enhancement) was 
first studied in the 1960s W Lax in the context of single¬ 
mode detuned gas lasers Q. The Lax a factor is 1a^, 
where ao is th® normalized detuning of the lasing fre¬ 
quency from the atomic resonance, i.e., ao = ; 

which is equal to the ratio of the real part of the gain 
permittivity to its imaginary part, or equivalently the ra- 

tio , where Arig is the refractive index change due 

to fluctuations in the gain. Two decades later, Henry 
derived an amplitude-phase coupling enhancement fac¬ 
tor of the same general type in semiconductor lasers [m, 

_ Re Artg latter case these refractive-index 

^ Im Arig ’ 

changes arise from carrier-density fluctuations and take 
a different form. Here, we are considering atomic gain 
media, so our a factor generalizes the Lax form. 

The difference between our single-mode generalized a 
factor dSOl) and that of Lax arises because we take into 
account spatial variation in the gain permittivity due to 
spatial hole-burning and also the non-Hermitian (com¬ 
plex) nature of the lasing mode. Hence we expect our 
factor to reduce to the Lax factor in some limits. For 


instance, consider the situation that was discussed in 
the last paragraph of Sec. HD of a low-loss Id Fabry- 
Perot cavity laser, operating near threshold. In this 
case, the nonlinear coupling coefficient is approximately 

(j ^ ^ show that the general¬ 
ized a factor is a = ^ ~ (the last approxima¬ 

tion is valid since in essentially all realistic cavities, the 
modes can be chosen to be predominantly real, i.e., have 
small imaginary parts). 

In many cases, however, our a deviates from the tradi¬ 
tional factor uq. An obvious example is when the lasing 
frequency precisely coincides with the atomic resonance 
frequency. In this case, the traditional factor vanishes, 
but a does not necessarily vanish. In the next section, 
we calculate and discuss the characteristic properties of 
the generalized a factor for two Id laser structures. 

B. Generalized single-mode a factor 

In this section, we evaluate the differences between the 
generalized and traditional a factors in Id model systems. 
We solve the full nonlinear SALT equations using our 
recent finite-difference frequency-domain (FDFD) SALT 
solver [ 2 l|. 

The generalized factor a can deviate significantly from 
the traditional factor ao when the latter is large (a similar 
argument was made in H). To see this, let us write 
the nonlinear coupling coefficient qualitatively as C oc 
(1 -I- tao)(l + 1/3), where the term 1 -|- iao is associated 

with the atomic lineshape-— and the term l-\-iB 

is a complex factor due to the remaining integral factors 
(we refer to the latter term as the modal contribution to 
the a factor). Typically /3 <C 1 and, consequently, the 
generalized factor is approximately a~ao + ,5(l + a§), 
so the difference between the generalized and traditional 
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factors grows quadratically with oq. 

To verify this argument, we study a model system in 
which the magnitude of ao can be controlled. Consider 
a quarter-wave dielectric photonic crystal (PhC), with 
a defect at the center of the structure (the geometry is 
depicted in the upper inset of Fig. [5^, similar to the struc¬ 
ture that was studied in Fig. S]). Adding enough layers 
of the periodic structure on each side of the defect to 
mimic an infinite structure, one finds that the system 
has a localized mode in the vicinity of the defect (lower 
inset), whose resonance freq uency is fixed to a real value 
within the energy gap [33|. To study finite-threshold 
lasers, we introduce gain and some passive loss (i.e., a 
positive imaginary permittivity term, which pushes the 
resonance poles away from the real axis in the complex 
plane). Since the resonance frequency of the defect mode 
is fixed by the geometry, by varying the resonance fre¬ 
quency of the gain, we control the detuning of the lasing 
mode from the atomic resonance, thus controlling the 
size of ao- As demonstrated in the figure, the deviation 
|a — ao| grows as the detuning Aiy = increases. 

The openness of the cavity also results in an enhance¬ 
ment of the a factor; the more open it is, the larger is 
the necessary imaginary part of the lasing mode, which 
causes a deviation from the standard formula. In order to 
test this prediction, we evaluate the generalized a factor 
for an open-cavity laser (Fig. [Sb), where we can control 
the radiative loss rate through the cavity walls and, con¬ 
sequently, this part of the modal contribution to a. We 
consider a cavity which consists of a dielectric slab (with 
permittivity £c) surrounded by air on both sides, with 
gain spread homogeneously inside the slab (upper right¬ 
most inset). The reflectivity of the cavity walls is deter¬ 
mined by the difference in cavity and air permittivities 
Ae = Ec — So- For relatively small dielectric mismatch, 
the cavity is relatively low-Q and our a factor differs sig¬ 
nificantly from the Lax factor. As Ae increases and the 
cavity Q increases, the generalized a factor converges to 
the original factor, so that the red and blue curves in the 
figure overlap. 

Unlike a photonic-crystal defect-mode cavity where 
there is a finite bandwidth of confinement [s^, this di¬ 
electric cavity has an infinite number of possible lasing 
resonances and thus when we sweep Ae, the a factor 
peaks periodically. This is because the free spectral range 
of the cavity is Auj k. and, therefore, changing 

ec corresponds to shifting the passive resonances and, 
consequently, the lasing modes. Every time a lasing 
mode crosses an atomic resonance, ao vanishes and corre¬ 
spondingly a becomes very small. The traditional factor 
is maximized when the atomic resonance is equidistant 
from two passive modes. The peak value is proportional 
to the free spectral range and, therefore, we find that it is 
proportional to 1/This type of effect may not have 
been observed previously because in macroscopic cavi¬ 
ties, the cavity resonances are very dense on the scale 
of the gain bandwidth, so the lasing mode can never be 
substantially detuned. However, in microcavities with 


large free spectral range, this could be an important ef¬ 
fect. Another intriguing property of the generalized a 
factor is that it varies discontinuously at the peaks (as 
is shown more clearly in the upper left-most inset). The 
traditional factor ao depends only on the mode detun¬ 
ing from resonance, so it approaches the same value on 
different sides of the peak. In contrast, the generalized 
factor a depends on the mode profile E^, which differs 
between the two interchanging laser modes on different 
sides of the peak, producing the observed asymmetry. 


C. Generalized multimode a factor 


Our multimode linewidth formula includes linewidth 
corrections from neighboring modes, which enter through 
the generalized a factor (since phase fluctuations in each 
of the modes couple to amplitude fluctuations in all 
other modes due to saturation of the gain). Accord¬ 
ing to the traditional ST formula ([2), when phase cross¬ 
correlations between different modes are neglected, each 
resonance-peak width is inversely proportional to the cor¬ 
responding modal output power. We find that when 
phase cross-correlations are included, the linewidth of 
each mode is a sum of inverse output powers of all the 
other modes—a type of multimode Schawlow-Townes re¬ 
lation. To see how this comes about, recall that the gen¬ 
eralized a factor, as given by (1551) . is proportional to 

1BA“^ R (1BA“^)^ . We show in appendix C that in¬ 


dividual factors in the product scale as [BA oc 
where Ojo is the steady-state amplitude of the j’th mode. 
Therefore, the multimode a factor is proportional to the 
2 (const) xijjj . 

sum: a)V) >,,■ -r-, i-e., a sum over terms which 

scale as inverse output powers. 

In the two-mode case, the linewidth formula for a las¬ 
ing mode in the presence of a neighboring mode is given 
explicitly by 


Zttio 


Pll 

■ /^I ■ 

^11^22 ^21^21 

' , R22 

■y^Ry^I /^I /^R ■ 

O 11 O 12 ^11^12 

2aio 

y^R y^R y^R y^R 

1^11^22 ^12^2lJ 

za2o 

y^R y^R y^R y^R 

1^11^22 ^12^2lJ 


(43) 


where = ReCi, and C?- = ImCi,. (A similar ex- 
pression was derived in [ 39 |, by using a phenomenologi¬ 
cal version of the two-mode TCMT equations.) As pre¬ 
dicted by the multimode ST relation, the last term in 
( 1431 ) is inversely proportional to the output power of the 
second mode o^q. This term becomes significant when 
the power in the first mode greatly exceeds the power in 
the second mode (i.e., when Pi ^ P 2 ), correcting the un¬ 
realistic Schawlow-Townes prediction that the linewidth 
vanishes when Pi —>■ 00 ; a similar argument was made 
in [ 3 ^. Fig. |2i presents the spectrum of a two-mode in¬ 
stantaneous model m in the parameter regime where 
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FIG. 6. (Color online) (a) Spectrum near a resonance peak in the presence of an additional mode. Numerical simulations of 
(I15II (red curve) and analytic single-mode (blue) and multimode (black) formulas. The simulation parameters are chosen so 
that there are two lasing modes with the same steady-state amplitudes ako — 1 and diffusion coefficients Rkk ~ 0.05, and with 
substantial cross correlations: Ckk = 5, Cki = 4 -|- 4i, fc yf f (in arbitrary frequency units), (b) Linewidth of central resonance 
peak vs. output power in the neighboring mode [aio = 1 and a 20 £ (0,3)]. Simulated spectrum (red) and analytic single-mode 
and multimode formulas (blue and black curves). The point aio = a 2 o = 1 is encircled, and corresponds to the parameter 
values of Fig. [6^. 


cross-correlations between the two modes are significant. 
The linewidth of the simulated spectrum (red curve) is 
in complete agreement with the generalized formula (1431) 
(black curve), but deviates substantially from the single¬ 
mode formula (1^ (blue curve). In order to reach the 
regime where this deviation is substantial, in practice, 
one needs to design a cavity in which the two lasing 
modes have comparable amplitudes and detunings from 
the atomic resonance frequency. 

Eq. (|43l) predicts an unphysical divergence near the 
second threshold, i.e., when 020 —f 0 (see black curve 
in Fig. EJd). In retrospect, this singularity is to be ex¬ 
pected, since the assumptions of our derivation break 
down in this limit. (No te that an equivalent diver¬ 
gence was present in @.) In calculating the phase 
variance, we assumed that amplitude fluctuations in all 
modes were small compared to the steady-state ampli¬ 
tudes (i5i <C Uio), and this assumption is no longer valid 
near threshold. The N-SALT TCMT equations (ITT]) are 
still valid, however—it is only their analytical solution 
for that is problematic. Therefore, we study the 

threshold regime numerically, via stochastic simulations 
of the N-SALT TCMT equations. As shown in Fig. [Bb, 
the simulated linewidth of the first mode approaches a 
finite value near the second threshold (red curve), and 
this value is significantly larger than the linewidth predic¬ 
tion one obtains when neglecting the second mode (blue 
curve). Even at the threshold, noise in the second mode 
mixes with the first mode through off-diagonal nonlinear 
coupling terms, thus increasing the linewidth. 

Linewidth enhancement at the thresholds of neighbor¬ 
ing lasing-modes suggests that the linewidth must also 
be enhanced below the modal thresholds [in the regime 


where radiation from non-lasing modes is incoherent, 
commonly called amplified spontaneous emission (ASF)]. 
We believe that this phenomenon could be explored using 
a future generalization of our formalism, with some mod¬ 
ifications (extending earlier work on linewidth 

enhancement from ASF). 


VII. FULL-VECTOR 3D EXAMPLE 

In order to illustrate the full generality of our approach, 
we apply it in this section to study a three-dimensional 
photonic-crystal (PhC) laser. The steady-state proper¬ 
ties of this system (i.e., the lasing threshold and mode 
characteristics) were previously explored in [2l|. We use 
those solutions here to calculate the laser linewidth [us¬ 
ing ([3|)] , and we compare the relative contributions of the 
various correction factors. 

The simulated PhC consists of a dielectric slab pat¬ 
terned by a hexagonal lattice of air holes (Fig. [7^). A 
defect is introduced by decreasing the radii of seven holes 
at the center of the structure [s^ , giving rise to a doubly- 
degenerate mode which is situated at the defect (spa¬ 
tially) and in the bandgap of the lattice (spectrally). We 
select the TE-like mode out of the degenerate pair by 
imposing even and odd reflection symmetry at a; = 0 and 
y = 0 respectively, as well as an even reflection symme¬ 
try at 2 = 0. Staying close to a potential experimental 
realization, we choose the pump profile to be uniform in¬ 
side the high-index dielectric near the defect region, and 
zero elsewhere. We solve the SALT equations using our 
scalable FDFD solver, and track the evolution of the first 
lasing mode upon increasing the pump strength from zero 
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FIG. 7. (Color online) Linewidth correction factors for a 3d 
PhC laser, (a) The PhC consists of a hexagonal lattice of 
air holes (with period a = 1 mm and radius 0.3 mm) in a 
dielectric medium with index n = = 3.4. The slab has a 

thickness of 0.5 mm, with air above and below, terminated 
by PML absorbers. A cavity is formed by seven holes of 
radius of 0.2 mm. The pump is non-zero in the hexagonal 
region for height 2mm in the z-direction. (Borrowed from 
(b) Schematics of a 4-level gain medium. Levels |1) 
and 12) form the lasing transition, with resonance frequency 
oJa = 1.5 mm“^ and polarization decay yx = 2.0 mm“^. The 
population decay rates are 701/712 = 723/712 = 10^ and 
712/7X = 10~^. The pump rate is varied in the range 
^/ 7 i 2 = 0.4... 2.2. (c) Traditional (dashed) and generalized 
(solid) correction factors, as defined in Table. H] The total cor¬ 
rection is defined as the product of the (traditional and gen¬ 
eralized) Petermann, a, bad-cavity, and incomplete-inversion 
factors. The x axis is the DpjDth, where Dp is the SALT 
effective pump parameter (see text) and Dth is the effective 
threshold pump. 


to five times the first-threshold value. 

Typically, realistic laser structures do not use 2-level 
gain media, but employ a more complex optical scheme 
which involves multiple levels and transitions in order to 
achieve significant inversion and depletion of the ground- 
state population. In this section, we apply our formalism 
to a 4-level gain medium (Fig.jTJo), using a generalization 
of SALT [^, which finds the stationary multimode las- 
ing properties of an A"-level gain medium. As shown 
in [73|, an A^-level system can be mapped into an effec¬ 
tive 2-level system, which obeys the (2-level) SALT equa¬ 
tions with renormalized pump (Dp) and atomic relax¬ 
ation rates (7||). Consequently, the linewidth of a 4-level 
laser will be given by our generalized formula ([3]) with 
the appropriately renormalized coefficients. By choosing 
the decay rate between the lasing transition levels (712 
in Fig. [7]3) to be much smaller than the decay rates into 


the upper (723) and out of the lower (701) states, we can 
achieve substantial inversion and ground-state depletion. 
Consequently, the incomplete-inversion factor is approx¬ 
imately risp Ri 1, close to typical measured values (47l |. 

Fig. [3: presents the traditional and new correction fac¬ 
tors (dashed and solid lines respectively), as defined in 
Table. H We find that those factors are relatively small 
for this system and, consequently, the deviations between 
the new and traditional factors are small. A small Peter¬ 
mann factor arises since the first lasing mode has a rela¬ 
tively high quality factor (i.e., the cold-cavity resonance 
pole is at wq = 1.725 — 0.00512z mm“^ with a quality 
factor of Q_~ 700, in agreement with experimental re¬ 
alization [ 83 )• Moreover, the cold-cavity resonance lies 
well within the gain bandwidth, resulting in small a and 
bad-cavity corrections. The generalized factor a (solid 
purple line) is obtained from from (I17I30I) . Deviations of 
a from the traditional factor ag = (dashed purple 

line) are due to modal contributions to the a factor (see 
Sec. |VI1B). The generalized Petermann factor (full blue 
curve) is compared against the traditional factor (dashed 
blue line), which is expressed in terms of the SALT mode 
(instead of the passive cavity mode). The cavity region 
is taken to be the entire high-index medium. (Note the 
the generalized and traditional factors agree at thresh¬ 
old). Both the Petermann and a factors increase the 
linewidth. However, the generalized and traditional bad- 
cavity factors (full and dashed red curves respectively) 
lead to linewidth reduction. 

Last, we evaluate the incomplete-inversion factor risp. 
The inversion D{x) is found from the SALT solutions of 
the effective 2-level system. The excited state population 
-^ 2 ( 2 } can be derived straightforwardly, using the results 
of [TSj as follows. Assuming that the populations in the 
non-lasing levels |0) and |3) are at steady-state, one can 
express those populations in terms of the populations in 
the lasing transition |1) and |2). Then, by invoking the 
density conservation condition, where n is 

the atom number density and Ni are the individual level 
populations, one hnds that the population in |2) is given 

by 


N2 


n + tD 
1 + r 


(44) 


where r = 1 -b -b Having obtained expres¬ 

sions for D and for N 2 ^ we have all that is needed to 
calculate the incomplete-inversion factor Usp. We de¬ 
fine the “linear incomplete-inversion factor” (rigp dashed 
green line) as the ratio , [i.e., both the excited- 

I—I 

state population (|4^ and the inversion are evaluated at 
D = Dp, neglecting hole-burning effects]. The “nonlin¬ 
ear incomplete-inversion factor” (usp solid green line) is 
defined in Table. |T| The nonlinear factor risp coincides 
with the linear factor Usp at threshold, but exceeds the 
traditional factor at higher pumps. We also plot the total 
linewidth correction, which is defined as the product of 
the (traditional and new) Petermann, a, bad-cavity, and 
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incomplete-inversion factors. 


VIII. CONCLUDING REMARKS 

We presented a generalized multimode linewidth for¬ 
mula, obtained from the N-SALT TCMT equations for 
the lasing mode amplitudes, which we derived start¬ 
ing from the Maxwell-Bloch equations and using the 
fluctuation-dissipation theorem to determine the statis¬ 
tical properties of the noise. Our generalized linewidth 
formula ([31) reduces to the traditional formula ([21) for low- 
loss cavities and simple lasing structures, but deviates 
significantly from the traditional theories for high-loss 
wavelength-scale laser cavities. By basing our derivation 
on the SALT steady-state lasing modes, it is possible 
to apply our formula to cavities of arbitrarily com plex 
geometry (e.g., photonic crystal or microdisk lasers Il5l - 
Il8l| l and arbitrary openness (e.g. random lasers |43|). 
Also, since SALT includes to high accuracy the effects 
of spatial hole-burning, our formula includes both gain 
saturation and the spatial variation of the gain permit¬ 
tivity well above threshold, plus all effects due to modal 
couplings. From a computational point of view it is im¬ 
portant to point out that our formula is analytical and 
can be evaluated immediately from the output of a nu¬ 
merical SALT calculation without any significant com¬ 
putational effort. A manuscript describing a brute-force 
numerical validation of our theory against numerical so¬ 
lution of the Maxwell-Bloch equations is currently being 
prepared . Given only the laser geometry, the pump¬ 
ing profile, and characteristic properties of the gain (i.e., 
its resonance frequency Wa and decay rate 7 _l), our for¬ 
mula enables linewidth calculation, including a general¬ 
ized a-factor and accounting for temperature variations, 
at a level of generality that was not possible before. This 
generality is most important, of course, in cases where the 
new result is substantially different than previous theo¬ 
ries, and it would be interesting to study laser cavities in 
which the discrepancy is as large as possible. 

One such case is that of lasers which contain excep¬ 
tional points (EPs) in their spectrum, which are points 
of degeneracy where two (or more) eigenfrequencies and 
eigenfunctions coalesce [gl, [s^. EPs in laser systems 
have been explored recently, both theoretically and 
experimentally [s^. At the EP, the modes become self- 
orthogonal and that causes the denominator of (|3]) to 
vanish and is already known to greatly enhance the Peter- 
mann factor [9^ . Since a similar denominator appears in 
the integrals defining our generalized a factor (I12I30|) , we 
expect that our a will differ substantially from previous 
results near an EP (and similarly for the inhomogeneous- 
temperature correction). 

An important and exciting addition to the theory 
would be a treatment of amplified spontaneous emis¬ 
sion (ASE) from inodes below threshold; we believe this 
can be achieved by deriving TCMT equations for below- 
threshold (passive) modes, in which there is no steady- 


state oscillation (generalizing previous ASE work which 
used simplified models [^ll^). Incorporating the ASE 
contribution to the spectrum will allow us to follow the 
noise through the lasing thresholds, correcting the un¬ 
physical divergence which was discussed in Sec. IVIIB. 
More importantly, treating below threshold ASE should 
allow an ab-initio theory of LEDs in arbitrary cavities 

Euture work could also incorporate several additional 
corrections that were not treated in this paper. Our 
derivation applies to isotropic materials described by a 
scalar permittivity e, but extension to anisotropic per¬ 
mittivity e, magnetic permeability (fl), and even bian- 
isotropic materials would be very straightforward (e.g., 
for an anisotropic £, the only change is that e'Fp' factors 
and similar are replaced by E • (eE) etcetera, as in [s^). 
As discussed in Sec. lYlC, we are also able to exploit 
our framework to analytically solve for the relaxation- 
oscillation side-peak spectra, and are currently preparing 
a manuscript presenting this analysis [s^. We believe it 
will be possible to extend our formalism to handle non- 
Lorentzian lineshapes arising from frequency dependence 
(correlations) in the noise within the laser linewidth (goI - 
[64 |. as also discussed in Sec. |lVl Instead of treating the 
noise spectrum S'F(a;) as a constant S'F(a;^), one needs to 
include a first-order correction , e.g., by Taylor expanding 
S-p{uj) around it might be convenient to fit S-p{ijj) to a 
Lorentzian matching the amplitude and slope at since 
the Eourier transform of a Lorentzian is an exponential 
that should be easy to integrate. Finally, as noted above, 
although our derivation was for the two-level Maxwell- 
Bloch equations, a similar approach should apply to more 
complex gain media (including multi-level atoms [^ . 
multiple lasing transitions, and gain diffusion [t^.) The 
N-SALT linewidth theory can be generalized to account 
for these laser models following along the lines of our 
approach here. 
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APPENDIX A: DERIVATION OF N-SALT TCMT 

In this appendix, we derive the TCMT equations for 
the lasing mode amplitudes. Our starting point is the 
Maxwell-Bloch equations , which describe the dy¬ 

namics of the electromagnetic field in a resonator inter¬ 
acting with a two-level gain medium: 

VxVxE-|-ecE = -47rP -|- Fs, (Al) 
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P = -i(wa - m)P - (A2) 

47r 

b = -7|| [Zip - D + 2^i(E • P* - E* • P)] , (A3) 

where E is the electromagnetic field, while P and D 
are the atomic polarization and population inversion. 
(From here on, for brevity, we refer to D as the “in¬ 
version.”) Wa is the atomic resonance frequency, and 
7 _l and 7 || are the population and inversion relaxation 
rates. Zip is the external pump, which determines the 
steady-state inversion, and Ec is the passive dielectric 
permittivity. The field, polarization and inversion are 
measured in their natural units: Cc = Pc = ^■\/7||7-l/( 2^) 
and dc = fi7_L/(47r5^) respectively, where g is the atomic 
dipole matrix element [l4l^ . We introduce sponta¬ 
neous emission noise by including a random source term 
Fs = 4Tr§ in written in the frequency domain as 

Fs(x, w) = —'i47ru;J(x, w), (A4) 

where J(x, w) is a random fluctuating current, and the 
correlations of Fs(x, lu) are given by the FDT. 

Steady-state ab-initio laser theory (SALT) handles the 
noise-free regime of the Maxwell-Bloch equations (i.e., 
Fs = 0) and reduces this set of coupled equations to a 
frequency-domain nonlinear generalized eigenvalue prob¬ 
lem for the electric field E (as reviewed in Sec. 1.1). 
When noise is introduced (Fs ^ 0), the cavity field is per¬ 
turbed from steady-state and the nonlinear permittivity 
is modified (Sec. 1.2). This gives rise to a restoring force 
(denoted Fnl), which we calculate in Sec. 1.3. The noise- 
driven field E is then found by integrating the Green’s 
function (derived in Sec. 1.4) over the noise terms Fs 
and Fnl- Finally, the TCMT equations are obtained by 
transforming back into the time domain (Sec. 1.5). 


1. Review of SALT 

We begin by reviewing the steady-state theory. In the 
SALT approach, the steady-state electromagnetic field is 
expressed as a superposition of a finite number of lasing 
modes: 

Eo(x, t) = J2 E^(x)a^oe-*^'‘‘, (A5) 

where Eo(x, t) denotes the steady-state field and a^o are 
the steady-state modal amplitudes. The lasing modes 
(x) are real frequency solutions of the nonlinear eigen¬ 
value problem 

[V X V X -u}leo{u}^,ao)] E^(x) = 0, (A6) 

with outgoing boundary conditions. The effective per¬ 
mittivity has a linear (passive) term Sc and a nonlinear 
(E-dependent) gain term: 

eo(w,ao)=ecH-— Do{ao). (A7) 

UJ - UJs, + 17_L 


The steady-state inversion Z7o(ao) [which is a notation 
shortcut for Z7o({E^}, {w^}, {a^o})] is given by 


-Do(ao) — 


D, 


1 + E 




M CJa)^+7 




(A8) 


To avoid possible confusion, note that in previous SALT 
works, the steady-state inversion was denoted by D and 
Dq was the external pump parameter, whereas in this 
work, Dq is the steady-state inversion and Dp is the ex¬ 
ternal pump parameter. 


2. Noise-driven Maxwell-Bloch equations 

In the presence of a small noise source, the electric 
field and polarization can be written as superpositions of 
the steady-state lasing modes with time-dependent am¬ 
plitudes a^(t) and bbt)'- 

E(x, t)=Y^ E^(x)a^(Z)e"*“’'^* 

P(x,t)=^P^(x)54t)e-“'^‘. (A9) 

Substituting the perturbation ansatz (IA9I) into the po¬ 
larization equation (IA2I) . we obtain 

(5^ -I- za;^fo^)P^ = -f(wa - ij±)bp,Pfj, - 

47r 

(AlO) 

Taking the Fourier transform and rearranging terms, we 
find 


^* Se^, (All) 

where we have introduced the shifted frequency uj = uj/j_ + 
fl and the Fourier-domain envelopes a^(n) = a^(w —w^), 
ZZ^(n) and Z7(f2). The asterisk * denotes a convolution. 

Next, consider Eq. (ED in the frequency domain 

V X V X E-(E-h47rP) = Fs. (A12) 

When the spacing between adjacent lasing modes is much 
larger than their linewidths, a noise source with fre¬ 
quency uj ke LOp, excites only the mode E^(x). Equiva¬ 
lently, the Green’s function can be approximated by the 
contribution of the single pole at (Note that we re¬ 
quire only that the peaks in the laser spectrum above 
threshold are non-overlapping; we do not require isolated 
resonances in the passive cavity spectrum.) Therefore, at 
frequencies w « we can substitute (lAllI) into (IA12I) 
and obtain an effective equation for the noise-driven field 
E^(x,w): 

[V X V X -u}'^e{uj,a)] E^(x,w) = Fs(x,a;), (A13) 
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where the effective permittivity e(w, a) is given by 


e(a;,a)E^(x,a;) = 


+ 


7-L 


-D * cin 


E^(x). 


UJ - UJa+ Z7_L 

(A14) 

The second variable of e{u},a) denotes the implicit de¬ 
pendence of e on the modal amplitudes through the 
Fourier transform of the inversion D. We calculate D 
explicitly in the next section. 


3. The atomic inversion 

The noise source Fs perturbs the modal amplitudes 
from steady state, causing a change in the atomic inver¬ 
sion D. We neglect dispersion corrections to D (which 
amounts to setting = 0 in (lAlOl) (s^j) as these correc¬ 
tions do not affect the linewidth formula to leading order 
in the noise [see discussion following (|9l) in the main text]. 
From (IA3I) and (lAlOl) . we obtain 


where e(a;,ao) is the steady-state SALT permittivity 
which was defined in (TATI) . Xi,{oJ,aQ) is the permittiv¬ 
ity differential due to deviation in the modal amplitude 
Qu [which we denote by C3)]' 

r. f 1^ 

— I ■ 1/ ^2| 2 I j / \ 

W-Wa + *7-L \{UJ^ - UJar + JlW 

(A20) 

and Aa^ is the Fourier transform of the time-averaged 
modal deviation from steady state 

Aa. = 7 (x) (A 21 ) 

Substituting the permittivity expansion (jA19l) into 
Maxwell’s equation (IA13I) . we obtain 

[V X V X -u}'^e{uj,ao)] E^(x,ca) = FNL(x,a;) + Fs(x,w), 

(A22) 

where the nonlinear restoring force is 


D = —f\\ D-Dp + j: 


7^ 




IE,. 


(A15) 


In order to solve (IA15I) . we linearize the time depen¬ 
dent products \a^\'^D in the sum around the steady state 
\a^\‘^D Ki alQDo + Do{\a^\'^-alg)+al{D-Do), where Dq 
is the steady state (SALT) inversion (lASp . To simplify 
the notation, we define the local decay rate 


7(x) = 7|| 



7i 


.)^+7i 



(A16) 

The second term in (IA16I) gives precisely the increased 
atomic decay rate due to stimulated emission. Using the 
definitions above , (IA15I) becomes 


D =-x(ic){D - Do) 

-7\\Do-Y1 


7i 


{UJ„ - LUa)^ + X'i 


|E.|2(|a.|2-|a.on, 

(A17) 


which we can integrate, and obtain 


D = Dq + D{ 


7i 


IE,. 


711 


{ujy - Wa)^ + 7i ' 
y* dt'e-'^(")(‘-‘')(|a.o|' - la.(t')l')- (A18) 


Having derived an explicit expression for D{t), we sub¬ 
stitute its Fourier transform D into the effective permit¬ 
tivity (IA14|) and obtain 


e(uj, 0 )% e(w, ao)E^ -I- ^ x,.(a;, ao)Aa^ * E^, (A19) 

V 


FNL(x,a;) =a;^^x,,(w,ao)Aai. *E^(x,w). (A23) 

V 

The left-hand side of (IA22I) is just the linearized steady- 
state equation (|A 6 | 1 . and the nonlinear correction to the 
effective permittivity due to the noise Fs appears as 
an additional source term Fnl- As noted above, the 
noise-driven field E^ is found by integrating the Green’s 
function of the steady-state equation (IA 6 |) over the noise 
terms Fml and F 5 . In the following section we derive an 
approximate formula for the Green’s function. 


4. The linearized steady state Green’s function 

The single-pole approximation of the Green’s function 
is valid for frequencies near the resonances a; Ri as long 
as the spectrum consists of non-overlapping resonance 
peaks, i.e., when the spacing between resonant modes 
exceeds the modal linewidths. First, let us rewrite the 
left-hand side of (IA6I) as an operator acting on the 
field E(x, w): 

£tjE(x, w) = (V X V X —uj^£o(‘^,ao)) E(x, w). (A24) 

Next, we choose a complete set (see below) of eigenfunc¬ 
tions E„(x, w) and eigenvalues A„(a;) of the operator £aj- 

£ujE„(x, uj) = An(uj)E„(x, uj). (A25) 

We define the inner product of two vector fields, A(x) 
and B(x), as {A,B) = f dx A(x) • B(x). The operator 
£aj is complex symmetric under this inner product, i.e., 
(A,£a,B) = (£ajA,B) [H, 113 . Therefore, we use un¬ 
conjugated inner products throughout the derivation. In 
order to treat the set {E„} as a discrete (countable) basis, 
a convenient theoretical trick is to place the system in a 
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box with absorbing boundary layers in which the absorp¬ 
tion turns on more and more gradually. This procedure 
also gives the states E„ finite norms {En,En). Because 
the operator is non-Hermitian, completeness of the ba¬ 
sis can break down at an “exceptional point” [^, [s^, 
but exceptional points are not generically present—they 
must be forced by careful tuning of parameters. There¬ 
fore, we assume completeness in this manuscript and will 
treat the influence of exceptional points (self-orthogonal 
modes) as a limiting case in a future paper, as discussed 
in Sec. IVIIIl 

Let G{uj, X, x') be the Green’s function of the operator 
Cu:, defined via CuiG{u},x,x') = 5 (x —x') [9l|. Given the 
complete set of eigenfunctions and eigenvalues iE„, Ap.j, 
the Green’s function can be expressed as the sum 

G(w, X, x') = ^^ 
n A„(w) • 




dx E^(x) 


(A26) 


Each lasing mode is associated with an eigenvalue 
of , which has a zero at a real frequency w = . Gon- 

sequently, G(uj, x, x') has a pole at oj^ and at frequencies 
near it is dominated by a single term in the sum. 
Expanding A^(aj) around the pole A^(a;) « (w — 

(where \'^= )) we obtain 


G^m(w,x,x') 


E;.(x)EJ(x') 


(w -a;^)A' 


dx E2(x) 


(A27) 


Green’s function over the source terms FNL(x^aJ) and 
Fs(x',a;): 

E^(x,a;) = ^a;^ / (ix'G(x, x', u;)xi,(u;^, ao)Aai, * E^ 
+ I dx'G{x,x',u})Fs. (A30) 


In the first term on the right-hand side, we approximate 


O' 


because the correction term is O 


{uj — uj^) ■ Aa 


which is second order in the noise. Substituting the 
single-pole approximation (jA29l) in (IA30I) yields 


Em^a* = 

1/ 



E^/x>.(w/.,ao)E2 ^ _ 

—p - - — Attiy * 

J da;E2(x) [uj%{uj)]'^ 


E^ /Fs(x',u;)E^ 


(A31) 


Finally, multiplying both side by w—and taking the in¬ 
verse Fourier transform, we arrive at the N-SALT TCMT 
equations, which govern the evolution of the modal am¬ 
plitudes a^: 



(A32) 


In order to evaluate A(^, let us rewrite Cuj as C^j ~ 
Ga;^-|-E(w), where = Vx Vx-a;^eo(w^) and E(w) = 
- [iv'^eoiu})]'^ ( oj — w^). According to the Hellmann- 
Feynman theorem, the derivative of the eigenvalue A^(a;) 
with respect to ui is given by 


A(, 


/ 


dx E2(x) [-uj%{u})]'^ 
Jdx E2(x) 


(A28) 


The nonlinear coupling coefficient is 


/(x) = —iu 


X^{Uf,,ao)El 


dxFlix) [u;2go(a;)]' 


and the Langevin force is 


U{t) = * 


. /[Fs(x',t)e-“A]E, 


dxFlix) [a;2eo(w)]' 


(A33) 


(A34) 


and substituting (IA28I) in (IA27|) , we find that for frequen¬ 
cies near the resonances a; ~ the Green’s function is 
approximately 


APPENDIX B: LINEWIDTH OF THE 
MULTIMODE TIME-DELAYED MODEL 




Ea.(x)EJ(x') 


(w^-w) / dxFlix) [a;¥o(w)]' 


(A29) 


In this section, we calculate the laser linewidth for the 
multimode time-delayed model by generalizing the solu¬ 
tion strategy of Sec.|V]in the text. We begin our analysis 
with the discretized time-delayed N-SALT TCMT equa¬ 
tion (1551) (repeated here for convenience): 


5. The N-SALT TCMT equations 

Having derived an expression for the Green’s function, 
the noise-driven field can be found by integrating the 


d 






Ik j 


dfl + /m- 


(Bl) 
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Following the approach of Sec. El A, we linearize (EUby 
expanding the mode amplitudes around their steady- 
state values: = (a^o + (where (5^ a^jo): and 

we omit the terms 0{6^). Then, we introduce additional 
variables 


e^ = lk f (B2) 

where is the time-averaged amplitude deviation of 
mode /i = 1... M from steady state at the spatial point 
k = 1... N. Having introduced the auxiliary variables 
the set of integro-differential equations (IBlIl turns 
into a linear system of ODEs, which we solve by applying 
several linear-algebraic transformations to obtain a com¬ 
pact expression for the covariance matrix, as described 
in detail below. 

Introducing the vector : the linear system 

of ODEs is conveniently written as 

<5^ = - ^(2a^oa.oRe[C'^J)C^ + /^, (B3) 

uk 

<i>^ = - ^(2a^oa.oIni[C^J)C^ + ^ (B4) 

yk 

= + (B5) 

To simplify the notation further, we introduce the MxM 
matrices and (fe = 1... N), with entries 

= 2a^oa.oRe[C^'J) (B6) 

=2a^oa.oIm[C'^J), (B7) 

and we rearrange the set of equations (IB3IIB5I) in a matrix 
form [compare with (I24I25P ]: 


where we have the introduced the shorthand notation 

^ E f ("')) (B13) 


for the second term, which is responsible for the general¬ 
ized a factor. 

In the remainder of this section, we calculate First, 
we solve the set of ODEs for and S (IB8IB10I ). and then 
we substitute the solution for S into (jBI3|) and evaluate 
the integrals. To this end, we begin by rewriting the 
equations for and S more compactly. We define the 
[(TV -I- 1) • M] X I vectors x and F and the [(TV -|-1) • TVf] x 
[(TV -I- 1) • TVf] matrix IK: 


(\) 


/fH\ 

e 

, F = 

0 



1 0 / 


K = 


/ 0 Ai A 2 ... Aat \ 
Ai -Ai 0 ... 0 

A 2 0 —A 2 0 

: : 0 0 

V Aat 0 ... 0 —Aat J 


(B14) 


where A^ are block-diagonal MxM matrices with jk on 
the diagonal entries, and the zeros in the definition of IK 
are block MxM zero matrices. Using these definitions, 
the equations for and 8 (IB8IB10I) can be conveniently 
written as 


—X. = —Kx -I- F. 
at 

The solution of (|B15p is 


(B15) 


= + (B8) 

= + (B9) 

= (BIO) 

The autocorrelation matrix of the phase vector #, 
which we calculate in this section, is determined by the 
autocorrelation matrix of the Langevin force 

(f(t)r^(t')) = 1R<5(T-T')- (BII) 

In order to compute we solve (IBOp by straight¬ 

forward integration. We find that the phase covariance 
matrix is a sum of a “pure” phase-diffusion term, pro¬ 
portional to and an amplitude-phase coupling term, 
proportional to J: 


,(t)= / 

J L J mp 

P 


(BI 6 ) 


and, in particular, the solution for is 




Mk-\-p,s 


mt')- (BI7) 


st 

For ease of notation, let us denote the {k + I) MxM 
block in the first column of the matrix ) by 

the shorthand notation [e 1 , so that $k = 

/p ^^'[ 6 “"^^*“* )]fe+i 4 f^(t'). Substituting the expression 
for into and using the autocorrelation function of 
the Langevin force (IBllI) . we obtain 


J 


1 




dt'dt"ds' 




fc-l-1,1 




IR 


[ 2 +'^) 1 ^ 1 ’ 

(B12) 

^ 2 

e ^ ^ 


J £+1,1 


(B18) 
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We proceed (not shown) by diagonalizing the matrix IK 
and evaluating the integrals in (IB18I) . (The intermedi¬ 
ate steps depend on the eigenvalues of IK and the matrix 
of eigenvectors, but the final result can be expressed in 
terms of the matrix inverse 1K“^). In the long-time limit, 
we keep the leading order term (which grows linearly in 
time) and we obtain 

(B19) 


In order to complete the derivation of the linewidth 
formula, we use the following identity: 


[K] 


-1 _ 
fc+1,1 “ 



(B20) 


which we prove below. Noting that 

Bfii, = and using the identity (IB20|) . we find that 

(IB19I) reduces to 

[R 'T 

J = BA"^ - (BA-i) , (B21) 


l^\ 

1 

VI/ 


-1 



Therefore, we obtain [IK]^.^^ 



-1 


APPENDIX C: LEMMA FROM SEC. VI.C 


[BA ^]ij OC 


ajQ 


In Sec. EDb, we present a multimode Schawlow- 
Townes relation, which states that the linewidths are 
proportional to a sum of inverse output powers of all 
the other modes. This result arises from a lemma which 
we prove here. We use the standard matrix-inverse for¬ 
mula [ 9 ^ 


where the adjugate matrix is defined as 

adjA= ((-1)*+^M,,E. (C2) 


which completes the derivation of the linewidth formula 
in the most general time-delayed model. In particular, 
and somewhat remarkably, the 7 terms completely can¬ 
cel in the computation of the first column of the matrix 
inverse, and drop out of the final result. 

Proof of the identity (IB20I1 : We use Schur comple- 
ment [92| for the lower-left corner of a matrix inverse: 

/ab\~E/ * 

V C D j “ -D-1C(A - BD-^C)-! * ) ’ 

(A and D need to be square matrices). Decomposing the 
matrix IK into the blocks 


M is the cofactor matrix, i.e. the matrix whose (f,j) 
entry is the determinant of the (i, j) minor of A (which 
is the matrix obtained from A by deleting the i’th row 
and the j’th column). From the definition of A (i.e., 
Aij = Re[Cij]aionjo), h follows that 


Mij — 



^iOO^joQij ■ 


(C3) 


where Qij (and later Q) denote constants that may de¬ 
pend on i and j, but are independent of the modal am¬ 
plitudes. Note also that 


A=(0) , B=(AiA 2 ... Atv ) , 



/ Al \ 


/ -Al 

0 ... 

0 \ 


A2 


0 

—A2 

0 

c = 


' D = 


0 ■■■ 

0 


Aat j 


1 0 

... 0 

—Kn / 


we can calculate the lower-left corner of [IK]j,^^ ^ 


det A = ]^a^o • Q. 

k 


Using (ICl]ICl) . we obtain 



Qij' 

O^iOOjjO 


(C4) 


(C5) 


Therefore, one can easily see that the lemma follows, 
since [BA-i],, oc a^ako ■ oc 


-D-iC(A- BD-^C)-^ = (D-iC)[B(D-iC)]-i = 






I 


( Al A 2 ... Aiv ) 

I 

llj 





APPENDIX D: COMPARISON WITH THE 
SCATTERING-MATRIX LINEWIDTH FORMULA 

In a recent scattering-matrix based linewidth the¬ 
ory 0 , Pillay et al. obtain a formula for the linewidth 
of a one-dimensional laser system, expressed in terms of 
integrals over the modes which solve the nonlinear SALT 
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equations. In this appendix, we prove that their for¬ 
mula (which applies to Id systems) is equivalent to our 
linewidth formula (jS]) (except that their formula gives a 
spatially averaged incomplete-inversion factor and omits 
the a factor). 

In the scattering-matrix approach, the lasing modes 
are described as purely outgoing wave functions ipo, 
which satisfy the nonlinear SALT equation 

V X V X tpoix) - ujQe{x,u}o)i’o{x) = 0, (Dl) 

and can be expressed as a superposition of outgoing chan¬ 
nel modes outside of the laser region 

tjjo{x) = '^bkUk{x,u!o) for r ^ C. (D2) 

k 

C denotes the scattering region (i.e., £ = 1 for a; ^ C). 
Note that 'ipo{x) is precisely the same the mode Eo(x) 
(which was used in Sec. m inside the cavity region.] The 
outgoing mode-amplitudes b are normalized to the value 
of ipo at the cavity boundary {x = L) 

h^h = ro{L). (D3) 


The apparent difference between our formula and the 
linewidth formula in is that the integral term in the 
denominator of our linewidth formula(l3l) is replaced by a 
sum of two terms in the scattering-matrix approach 



space 


ewo 


ujq de 
2 c/cJn 


V'o 


ih^ b 


+ dx 

Jc 


Wn de 


V'o- 


(D4) 


In order for the two formulas to agree, we need to show 
that 


r°° ib^b 

Wo J dxipl^x) = (D5) 

(where we have used the fact that £ = 1 outside the 
cavity region). We show that the latter condition (ID5I) 
holds for any solution tpo of jnn) which satisfies outgo¬ 
ing boundary conditions. One way to impose outgoing 
boundary conditions is to invoke the limiting-absorption 
principle (i.e., add loss to eliminate incoming waves from 
infinity and take the limit of infinitesimal absorption at 
the end of the calculation [allll ). Formally, we define 
the integral on the left-hand side of (IDSI) as 

poo nOO 

/ dx4’1{x) = lim / dxe~^^il)1{x). (D6) 

Jl s-»-o+ Jl 


By substituting V'o(a^) = into (|D6I) and taking the 
limit of s —>■ O'*", we obtain J^dxijjoix) = = 

Ajdj^b, and since ujq = ckg this finishes the proof of 
(ID3|) (with the units convention of c = 1). 


APPENDIX E: ZERO-POINT FLUCTUATION 
CANCELLATION 


The hyperbolic cotangent factor in the FDT dT51) arises 
as a sum of a Bose-Einstein distribution and a 1/2 fac¬ 
tor stemming from quantum zero-point (ZP) fluctua¬ 
tions [zainl, and this is why it does not vanish in the 
limit of zero temperature (/3 —>■ 00). However, it turns out 
that contribution of this ZP term cancels in the linewidth 
formula, as was shown by Henry and Kazarinov from 
a quantum-operator viewpoint, and it is convenient to ex¬ 
plicitly subtract the ZP term from the hyperbolic cotan¬ 
gent as in ([S]) and (|2T]| . Here, we provide a purely classi¬ 
cal explanation for why this cancellation occurs, and why 
it is important to perform the explicit subtraction in or¬ 
der to eliminate a subtlety arising from the definition of 
outgoing boundary conditions. 

The FDT has a hyperbolic cotangent factor, and when 
we apply the FDT to find the correlation function 

in Sec IV, the same hyperbolic cotangent factor arises in 
the R integral, appearing in the form 


J da;|E^pIm£o(w^) • icoth 

= y da;|E^pIm£o(w^) Q coth 





for a lasing mode /r, where we have trivially added 
and subtracted the ZP 1/2 factor from coth. Now, we 
wish to analyze the final 1/2 term, which is the integral 
i / da;|E^pIm £o(a;;_j). Before we treat outgoing bound¬ 
ary conditions, let us consider the simpler case of a laser 
surrounded by an explicit absorbing medium, as in [^ . 
(This is also the situation in more recent computational 
models, for which one uses a finite spatial domain sur¬ 
rounded by absorbing layers [2l|.) For any steady-state 
lasing mode (real w^), the net gain -|- loss is zero, but 
da;|E^pIm £0(0;^^) is proportional to the net power 
absorbed or gained by the electric field [9^ and hence 
this integral is zero. Therefore, in such a case, whether 
or not we include the 1/2 factor is irrelevant, because the 
±1/2 terms integrate to zero. 

However, a subtlety arises in this integral in the com¬ 
mon case where the laser is surrounded by an infinite 
zero-temperature (/3 = 00) lossless medium with outgo¬ 
ing radiation boundary conditions. Outgoing boundary 
conditions can be defined mathematically by the limit¬ 
ing absorption principle 113 ■ takes the lossless 
medium to be the limit of a lossy medium as the losses 
go to zero from above, which can be expressed by writ¬ 
ing £ as £ -I- . Just as in appendix D, the correct 

approach is to take the lossless limit after solving the 
problem, i.e. the 0“*" limit is taken outside of the inte¬ 
gral. Before we take this limit, it makes no difference 
whether the 1/2 factor is included, just as above: it inte¬ 
grates to zero. However, after we take the lossless limit, 
there is no explicit absorbing region (Im £ > 0) in the 
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integral (the absorption has been “moved to infinity” in 
some sense), so if we perform the coth integral without 
subtracting 1/2 then we would obtain an incorrect con¬ 
tribution from the ZP fluctuations in the gain medium 
(which should have been canceled). Instead, if we inte¬ 
grate against 1/2 coth —1/2, the result is correct without 
requiring any explicit contribution from the absorbing 
boundary conditions. 

Note that if the laser is surrounded by an infinite loss¬ 
less medium at a positive temperature, then there is a 
nonzero contribution of incoming thermal radiation to 


the linewidth This can be included in one of two 

ways. In practice, we typical solve the SALT equations 
in a finite computational box with an explicit absorbing 
region, in which case no modification to our linewidth 
formula is required: one simply assigns the ambient tem¬ 
perature to the absorbing region. If, on the other hand, 
the outgoing boundary conditions are imposed in some 
other wa y ( e.g. semi-analytically as in earlier SALT 
work |19l . I20II4311971 ), then an explicit source term must 
be added to account for incoming thermal radiation, as 
in previous works [o^. 
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